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Introduction 


Energy-  and  intensity-modulated  radiotherapy  (MERT)  is  the  optimization  of  many  small  beamlets  of 
different  energies  and  intensities  to  produce  a  highly  conformal  dose  distribution.  Each  individual 
beamlef  s  dose  distribution  must  be  calculated  in  a  computed  tomography  based  representation  of  the 
patient.  The  relative  weight  of  each  beamlet  is  optimized  using  a  computer  implemented  mathematical 
technique.  When  the  optimized  beamlet  weights  are  applied  to  the  patient,  a  conformal  dose 
distribution  is  achieved.  This  is  called  inverse  planning  and  produces  very  conformal  dose 
distributions  when  coupled  with  MERT.  When  applied  to  radiotherapy  breast  treatment,  MERT  dose 
distributions  conform  to  the  treated  breast  while,  at  the  same  time,  spares  the  lung,  heart,  and 
contralateral  breast  to  excessive  irradiation.  Two  effects  that  thwart  MERT  treatments  are,  inaccurate 
dose  calculation,  and  the  effects  of  patient  organ  motion  and  daily  setup  uncertainty.  The  overall  aim 
of  this  project  is  to  use  accurate  Monte  Carlo  dose  calculation  as  the  basis  for  the  development  and 
implementation  of  mathematical  models  of  both  organ  motion  and  setup  uncertainty.  Monte  Carlo  is 
generally  accepted  by  the  medical  physics  community  as  the  most  accurate  dose  calculation  model 
available.  Once  Monte  Carlo  is  combined  with  realistic  models  of  organ  motion  and  setup  uncertainty, 
we  can  include  those  effects  into  the  calculation  of  beamlets  for  optimization.  Models  of  organ  motion 
and  setup  uncertainty  are  developed,  implemented,  and  experimentally  tested,  so  we  can  study  the 
improvement  in  the  MERT  dose  distributions  when  these  deleterious  effects  associated  with 
radiotherapy  treatment  are  included  in  the  beamlet  calculation  and  optimization  processes.  The 
successful  completion  of  this  project  will  facilitate  the  implementation  of  MERT  as  a  realistic  clinical 
option  for  breast  cancer  treatment.  All  patient-related  factors  should  be  included  in  the  dose 
calculation  and  optimization  processes  to  obtain  clinically  realistic  dose  uniformity  in  the  target 
volume  and  normal  tissue  sparing.  We  believe  that  that  the  inclusion  of  organ  motion  and  setup 
uncertainty  into  accurate  Monte  Carlo  dose  calculations  will  significantly  improve  optimized  MERT 
treatment  plans.  The  degree  of  improvement  in  the  final  treatment  plan  when  organ  motion  and  setup 
uncertainty  are  included  in  the  dose  calculation  of  beamlets  for  the  optimization  process  can  now  be 
investigated  after  the  successful  completion  of  our  specific  aims  for  the  first  year  of  research. 
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Body 

Traditional  radiotherapy  treatment  for  breast  cancer  consists  of  two  tangential  photon  beams.  The 
tangential  beam  configuration  is  used  to  limit  the  dose  to  the  lung  (and  heart  for  left  breast  treatment) 
by  limiting  the  beam  divergence  into  these  normal  tissues.  However,  even  with  this  approach, 
conventional  radiotherapy  breast  treatment  is  associated  with  unwanted  radiation  complications1'6. 
These  include:  high  dose  to  the  heart  and  lung  (resulting  in  radiation  pericarditis  and  radiation 
pneumonitis),  increased  radiation  scatter  to  the  untreated  breast  (resulting  in  an  increased  risk  of 
secondary  cancer),  and  poor  cosmetic  results  for  patients  with  larger  breasts.  A  novel  approach  using 
electron  beams,  rather  than  photon  beams,  has  been  proposed  to  overcome  these  problems7.  In  this 
technique,  a  conventional  broad  electron  beam  (25x25  cm2)  is  subdivided  into  many  small  electron 
beams  (lxl  cm2),  called  beamlets.  Each  beamlet  can  then  be  subject  to  energy-  and  intensity 
modulation.  This  is  done  by  initially  calculating  the  dose  distribution  for  each  beamlet  of  different 
electron  beam  energies  in  the  patient.  Then,  a  computer  implemented  optimization  algorithm  is 
applied  to  each  beamlet  in  order  to  determine  the  optimal  intensity  for  each  beamlet  energy.  The  result 
is  a  highly  conformal  dose  distribution  to  the  treated  breast  with  minimal  dose  to  adjacent  normal 
tissues.  This  approach  is  called  energy-  and  intensity  modulated  radiotherapy  (MERT)  and  is  a 
technique  called  inverse  planning.  It  is  the  physical  characteristics  of  electron  beams  that  give  MERT 
a  real  advantage  over  photon  beams  in  treating  breast  cancer  with  radiation.  Until  recently,  various 
theoretical  and  technical  limitations  have  not  allowed  the  clinical  fruition  the  dosimetrically 
advantageous  MERT  treatments.  Some  of  these  are  currently  being  overcome  with  other  DOD  funded 
research8,9.  Additional  theoretical  limitations,  which  are  the  subject  of  this  research,  is  the  organ 
motion  and  daily  setup  uncertainty  associated  with  radiotherapy  breast  cancer  treatment.  By  including 
the  effects  of  organ  motion  and  setup  uncertainty  into  accurate  Monte  Carlo  beamlet  calculations  for 
the  optimization  algorithm,  we  will  ensure  that  beamlets  calculated  for  MERT  will  be  optimized  for  a 
clinically  real  situation  and  thus,  the  planned  and  delivered  doses  will  be  the  same.  Current  techniques 
to  account  for  organ  motion  and  setup  uncertainty  is  to  simply  ignore  them  altogether  or  to  add  a  large 
margin  around  the  treated  volume.  An  increase  in  margin,  however,  results  in  an  increase  in  the 
volume  of  normal  tissue  irradiated,  which  in  turn,  increases  normal  tissue  complications.  The  MERT 
technique  will  be  a  very  promising  modality  for  breast  cancer  treatment  provided  all  clinical 
circumstances  are  accounted  for  in  the  planning/dose  calculation  process.  The  successful  completion  of 
the  project  will  help  the  clinical  realization  of  MERT  to  improve  the  target  dose  uniformity  and 


5 


A 


conformity  and  reduce  dose  to  the  lung,  heart  and  contralateral  breast. 

This  project  has  4  specific  aims:  (1)  develop  theoretical  models  of  organ  motion  and  setup  uncertainty, 
(2)  implement  the  models  into  Monte  Carlo  dose  calculations,  (3)  verification  of  the  new  Monte  Carlo 
programs,  and  (4)  evaluate  MERT  treatment  plans  for  10  breast  cancer  cases.  We  report  on  the 
research  accomplishments  associated  with  the  tasks  outlined  in  the  approved  “Statement  of  Work” 
below  for  our  research  between  August  1 , 2000  and  July  3 1  2001 . 

1.  The  Monte  Carlo  dose  calculation  code  system 

We  have  developed  a  general  purpose  Monte  Carlo  EGS410  user  code  MCDOSE1 1  that  can  be  used  for 
patient  specific  electron  dose  calculations.  Excellent  agreement  was  achieved  between  the  MCDOSE 
results  and  measurements12  and  the  well  benchmarked  EGS4/BEAM/DOSXYZ  code13.  Our  Monte 
Carlo  code  system  includes  the  option  to  produce  beamlets  for  both  photon  and  electron  treatment 
optimization.  The  requirement  of  dose  calculation  accuracy  for  truly  optimized  treatment  plans  has 
been  demonstrated  and  is  shown  in  Figure  1.  Calculation  inaccuracies  have  also  been  demonstrated 
with  inverse  planning  for  photon  beams14.  The  code  system  has  also  been  developed  and  demonstrated 
to  produce  beamlet  dose  calculations  for  MERT  optimized  breast  treatment  plans15. 

2.  Computer  model  of  organ  motion 

Tumor  motion  due  to  respiration  during  radiation  treatment  will  lead  to  “motion  artifacts”  in  dose 
distributions.  Intra-fraction  organ  motion  creates  a  broadened  dose  penumbra,  which  limits  the 
conformal  dose  coverage  to  the  tumor.  We  have  developed  analytical  models  to  trace  the  motion  of  the 
tumor  in  lung  based  on  the  initial  patient  computed  tomography  data  and  the  motion  data  measured 
during  a  radiotherapy  breast  treatment.  This  allows  the  reconstruction  of  realistic  dose  distributions 
received  by  the  patient  during  a  fraction  of  treatment.  The  patient  geometry  at  any  time  during  a 
breathing  cycle  is  reconstructed  using  the  computed  tomography  data  taken  at  the  inspiration  and 
expiration  phases  and  the  chest  wall  motion  measured  using  an  optical  motion  detector.  Accurate 
correlation  between  the  voxels  in  the  inspiration  (or  expiration)  geometry  and  the  voxels  in  the 
reconstructed  geometry  at  any  point  of  the  breathing  cycle  is  established  so  that  the  dose  to  a  voxel  (or 
an  organ)  can  be  accumulated  accurately  during  a  treatment16,17. 
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3.  Computer  model  of  setup  uncertainty 


The  intra-fraction  setup  uncertainty  is  a  random  phenomena1819  having  a  Gaussian  distribution.  The 
variance  of  the  Gaussian  distribution  is  different  for  each  translational  axis  (see  Figure  2)  and  depends 
on  the  anatomical  region  being  treated.  We  used  a  “real-time”  method  to  incorporate  setup  uncertainty 
into  our  Monte  Carlo  calculations.  This  was  done  by  perturbing  the  incident  particles  before  they  enter 
the  patient.  The  perturbation  was  done  by  sampling  from  three  individual  Gaussian  distributions,  one 
for  each  translational  x,  y  or  z  axis  in  the  patient  (i.e.,  Anterior-Posterior,  Right-Left  and  Superior- 
Inferior).  The  variance  of  the  each  Gaussian  distribution  was  obtained  from  the  literature19  and  is  in 
the  process  of  being  verified  by  measurement  in  our  department.  Setup  uncetainty  was  assumed  to 
consist  of  x,  y  and  z  translations  that  are  independent.  The  model  is  shown  schematically  in  Figure  3. 
Our  approach  to  include  setup  uncertainty  is  different  than  a  simple  Gaussian  smoothing  after  the  dose 
has  been  calcualated.  The  “real-time”  perturbation  of  each  incident  particle  before  it  enters  the  patient 
ensure  that  each  particle  that  enters  the  patient  “sees”  the  correct  source-to-surface  distance  and 
heterogeneities  within  the  patient.  For  a  Gaussian  smooting,  after  the  dose  calculation  is  completed, 
dose  near  heterogeneities  is  considered  in  only  an  approximate  way.  Ongoing  work  for  this  project  is 
to  compare  our  method  with  the  Gaussian  smoothing  method  of  the  static  dose  distribution  for  any 
differences  in  the  results. 

4.  Implementation  and  effects  of  organ  motion 

The  method  has  been  implemented  for  Monte  Carlo  treatment  planning  by  interfacing  with  the 
MCDOSE  code.  MCDOSE  phantom  files  for  different  phases  of  the  breathing  cycle  were  generated 
and  the  3D  dose  data  were  obtained  from  MCDOSE  simulations.  The  final  dose  distribution  was 
reconstructed  from  these  3D  dose  data.  The  purpose  of  this  work  was  to  compare  the  effects  of  organ 
motion,  specifically  thoracic  movement  on  conventional  photon  and  modulated  electron  beam  dose 
distribution.  The  Monte  Carlo  MCDOSE  treatment  planning  system  was  used  for  this  study  to 
investigate  breast  treatment  using  parallel  opposed  6MV  tangential  beams  and  MERT.  The  MERT 
treatment  consisted  of  applying  multiple  beams  consisting  of  6,  9,  12,  16,  20  MeV  energies  to  conform 
the  prescription  dose  to  the  target.  For  the  patient  computed  tomography  scans  used  for  initial 
treatment  planning,  the  dose  was  normalized  to  50  Gy.  Two  more  computed  tomography  slices  were 
generated  simulating  patient  exhalation  and  inhalation.  Treatment  planning  was  done  on  all  three 
images  using  the  same  parameters. 
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There  was  not  a  significant  change  in  dose  to  the  target  with  simulated  motion  of  one  cm  posteriorly 
and  anteriorly  with  respect  to  the  original  patient  computed  tomography  scan.  MERT  showed  an 
improvement  in  conformity  of  overall  dose  to  the  target  as  well  as  avoiding  dose  to  the  lung.  For 
photon  beam  treatment  the  right  lung  beneath  the  target  breast  shows  a  disparity  in  absolute  dose  and 
the  volume  of  tissue  irradiated  at  different  phases  of  the  breathing  cycle.  During  exhalation  7%  of  the 
right  lung  volume  received  a  dose  of  17.5  Gy  while  at  inhalation  7%  of  the  right  lung  volume  received 
a  dose  of  35.5  Gy.  At  mid-plane  between  the  two  extremes  of  the  cycle  only  2%  of  the  right  lung 
received  a  dose  of  35.5  Gy.  For  MERT,  the  lung  dose  remained  negligible,  however,  there  was  an 
increase  in  lung  dose  during  inhalation.  Ten  percent  of  the  lung  received  a  dose  of  5  Gy,  which  was  a 
dose  increase  of  5%  from  the  median  plane  of  breathing  cycle.  The  maximum  dose  to  the  lung  during 
inspiration  was  87%  of  the  prescription  dose  while  at  median  phase  of  the  breathing  cycle  the 
maximum  dose  decreased  to  80%  of  the  prescription  dose  and  at  expiration  this  dose  was  less  than  6% 
of  the  prescribed  dose  when  using  photon  beams. 


5.  Implementation  and  effects  of  setup  uncertainty 

Results  obtained  using  our  approach  to  include  setup  uncertainty  in  Monte  Carlo  IMRT  dose 
calculations  was  demonstrated  for  a  clinical  photon  beam  inverse  planning  case.  We  present  a  case  for 
the  treatment  of  a  metastatic  lesion  at  T1 1  -  T12.  This  patient  was  planned  with  9  coplanar  intensity 
modulated  fields  of  15  MV  photons.  Figure  4  shows  the  axial  dose  distribution  for  Monte  Carlo 
calculations  with  and  without  setup  uncertainty  for  the  first  clinical  case.  In  general,  the  dose  results 
look  similar.  However,  there  is  a  subtle  difference  in  target  coverage  and  slightly  larger  difference  in 
spinal  cord  dose.  The  difference  in  the  two  dose  distributions  (with  and  without  setup  uncertainty)  is 
made  clear  by  plotting  the  cummulative  dose  volumes  histograms  as  shown  in  Figure  5.  There  is 
similar  target  coverage  but  an  increased  dose  was  delivered  to  adjacent  tissues  including  the  spinal 
cord.  Other  critical  structures  are  distant  to  the  target  and  the  DVHs  with  and  without  setup 
uncertainty  are  deemed  to  be  clinically  the  same.  It  is  also  worth  noting  that  the  DVH  to  the  entire 
patient’s  body  (skin)  are  the  same  for  both  cases  so  the  total  energy  deposited  in  the  patient  is  the  same 
for  our  simulations  with  and  without  setup  uncertainty.  The  dose  is  simply  spread  over  a  large  volume 
to  other  tissues. 
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A  second  case  was  treated  for  a  lesion  in  the  thorax  and  planned  with  5  inverse  planned  coplanar 
fields  of  4  MV  photons.  An  interesting  aspect  is  demonstrated  in  the  case  shown  in  Figures  6  and  7. 
Figure  6  shows  an  axial  dose  distribution  for  Monte  Carlo  calculations  with  and  with  out  setup 
uncertainty  for  the  second  clinical  case.  In  Figure  7,  the  target  DVHs  for  the  second  clinical  case  show 
a  slight  increase  in  target  dose  homogeneity  when  setup  uncertainty  is  considered  in  the  IMRT  dose 
calculation.  We  believe  the  reason  for  this  is  that  the  hot-spots  within  the  target  will  be  smoothed 
away.  It  is  worth  noting  that  IMRT  treatment  plans  typically  have  a  large  dose  gradient  across  the 
target  volume.  The  gradient  can  be  upwards  of  30%  in  some  cases.  Setup  uncertainty  then,  will 
smooth  these  high  dose  regions  and  provide  a  more  homogeneous  dose  to  the  target.  However,  this 
smoothing  will  come  at  the  expense  of  an  increased  dose  to  adjacent  critical  structures.  The  amount  of 
that  dose  increase  may  depend  on  the  shape  and  location  of  the  target.  For  critical  structures  distal  to 
the  target  and  in  low  dose  regions,  the  effect  of  setup  uncertainty  is  not  large.  The  extent  to  which 
these  effects  are  reduced  when  using  MERT  is  being  determined. 

Similar  to  a  simple  Gaussian  smoothing  of  the  static  dose  distribution,  our  method  is  equivalent  to 
including  the  effects  of  setup  uncertainty  over  an  infinite  number  of  small  fractions20"22.  Our  method 
of  including  the  setup  uncertainty  implies  that  all  the  statistical  fluctuations  associated  with  setup 
uncertainty  have  been  averaged  out  to  the  maximum  possible  extent.  There  may  be  cases  with  less 
differences  and  there  may  be  cases  with  greater  differences  in  the  DVHs.  We  feel  that  by  including 
setup  uncertainty  in  our  method  at  least  gives  the  clinician  a  reasonable  estimate  of  the  possible  effects 
of  setup  uncertainty  during  treatment.  Additionally,  our  method  does  combine  proper  consideration  of 
a  patient’s  heterogeneities  and  irregular  contours  by  accurate  Monte  Carlo  dose  calculation  . 
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Key  Research  Accomplishments 


•  Develop  computer  model  for  organ  motion :  A  model  of  organ  motion  has  been  developed.  This 
model  uses  clinical  computed  tomography  scan  of  actual  breast  cancer  patients.  Thus,  the  models 
are  realistic  estimates  of  patient  positioning  and  immobilization. 

•  Implement  organ  motion  model  into  Monte  Carlo  code'.  We  have  implemented  the  organ  motion 
model  into  our  Monte  Carlo  code  called  MCDOSE.  This  model  is  implemented  in  a  general  sense 
and  can  also  be  applied  to  conventional  treatments  as  well  as  MERT  treatments  that  are  calculated 
with  MCDOSE. 

•  Develop  computer  model  of  setup  uncertainty'.  We  have  developed  a  model  of  setup  uncertainty 
based  on  published  data.  Since  it  has  been  shown  that  setup  uncertainty  is  a  random  process  in 
radiotherapy,  our  model  has  been  implemented  based  on  the  Gaussian  distribution.  The 
implementation  allows  the  adjustment  of  the  Gaussian  parameters  as  the  setup  uncertainty 
measurements  become  more  accurate  in  the  future. 

•  Implement  setup  uncertainty  model  into  Monte  Carlo  code :  We  have  implemented  our  model  of 
setup  uncertainty  in  our  MCDOSE  Monte  Carlo  code.  It  has  been  implemented  in  such  a  way  that 
it  can  be  “tumed-on”  and  “turned-off’  for  successive  Monte  Carlo  calculations  to  investigate  its 
effects  on  dose  results  and  optimization. 

•  Initial  study  of  the  effects  of  organ  motion  and  setup  uncertainty  on  optimization :  Our  initial  results 
show  that  the  effects  of  setup  uncertainty  limit  the  degree  of  optimization  of  the  dose  distribution  if 
setup  uncertainty  is  not  included  in  the  beamlet  simulation  for  optimization. 
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Reportable  Outcomes 

Peer-reviewed  papers  resulting  from  or  supported  in  part  by  this  grant : 

•  C.M.  Ma,  T.  Pawlicki,  M.C.  Lee,  S.B.  Jiang,  J.S.  Li,  J.  Deng,  E.  Mok,  B.  Yi,  G.  Luxton  and  A.L. 
Boyer,  Energy-  and  intensity-modulated  electron  beams  for  radiotherapy,  Phys.  Med.  Biol.  (2000) 
45:2293-2311 

•  J.S.  Li,  T.  Pawlicki,  J.  Deng,  S.B.  Jiang,  E.  Mok  and  C.M.  Ma,  Validation  of  a  Monte  Carlo  dose 
calculation  tool  for  radiotherapy  treatment  planning,  Phys.  Med.  Biol.  (2000)  45:  2969-2985 

•  Pawlicki  T,  Jiang  SB,  Li  JS,  Deng  J,  Ma  CM,  Monte  Carlo  dose  calculation  and  setup  uncertainty, 
To  be  submitted  to  Med  Phys,  May  2002. 

Meeting  abstracts  resulting  from  or  supported  in  part  by  this  grant : 

•  Ding  M,  Deng  J,  Li  JS,  Lee  MC,  Jolly  J,  Pawlicki  T  and  Ma  C-M,  Investigation  of  dose  correlation 
for  organ  motion,  Accepted  for  oral  presentation  at  the  43rd  Annual  Meeting  of  the  American 
Association  of  Physicists  in  Medicine  in  Salt  Lake  City,  Utah  from  July  22-26, 2001. 

•  Jolly  J,  Ding  M,  Pawlicki  T,  Chu  J  and  Ma  C-M,  Monte  Carlo  simulation  of  the  effect  of  thoracic 
motion  on  breast  treatment,  Accepted  for  oral  presentation  at  the  43rd  Annual  Meeting  of  the 
American  Association  of  Physicists  in  Medicine  in  Salt  Lake  City,  Utah  from  July  22-26, 2001 . 

•  Pawlicki  T,  T.M.  Guerrero,  S.B.  Jiang,  J.  Deng,  J.  Li,  A.L.  Boyer  and  C.-M.  Ma,  The  role  of  set-up 
uncertainty  in  intensity  modulated  treatment  planning.  Oral  presetation  at  the  42nd  Annual  Meeting 
of  the  American  Society  for  Theraputic  Radiology’  and  Oncology,  Boston,  MA,  October  22  - 
October  26, 2000. 

•  Pawlicki  T,  Jiang  SB,  Li  JS,  Deng  J  and  Ma  C-M,  Including  setup  uncertainty  in  Monte  Carlo  Dose 
Calculation  for  IMRT,  Oral  presentation  at  the  2000  World  Congress  on  Medical  Physics  and 
Biomedical  Engineering  (in  conjunction  with  the  42nd  Annual  Meeting  of  the  American  Association 
of  Physicists  in  Medicine),  Chicago,  IL,  July  23-28,  2000. 

Funding  applied  for  based  on  work  resulting  from  or  supported  in  part  by  this  grant: 

None 
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Conclusions 


v 


We  have  made  significant  progress  during  our  first-year  investigation.  We  have  successfully 
performed  the  tasks  scheduled  in  the  “Statement  of  Work”.  We  have  developed  models  of  organ 
motion  and  setup  uncertainty  for  implementation  into  our  MCDOSE  Monte  Carlo  code  and  these 
models  have  been  successfully  implemented  into  that  code.  Preliminary  results  indicate  that  if  organ 
motion  and  setup  uncertainty  are  not  included  in  the  pre-optimization  beamlet  simulations,  then  the 
beamlet  weights  and  subsequent  dose  distributions  are  not  truly  optimized.  Further  studies  as 
scheduled  in  the  “Statement  of  Work”  are  needed  to  complete  the  verification  of  our  implementation 
and  also  to  determine  the  degree  of  efficacy  of  MERT  treatments  for  breast  cancer.  In  the  process  of 
this  work  we  have  found  that  our  models  and  implementation  can  be  utilized  to  improve  the 
conventional  treatment  of  breast  cancer  as  well. 
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Appendices 

List  of  Figures  quoted  in  the  body  of  text: 

Figure  1 :  These  are  two  plans,  one  optimized  based  on  FSPB  beamlets  and  one  based  on  Monte  Carlo 
beamlets  (a).  All  other  parameters  are  held  constant.  The  result  is  more  conformal  and  homogeneous 
dose  distribution  to  the  target.  The  DVHs  show  better  coverage  and  critical  structure  sparing  (b).  The 
solid  line  denotes  the  Monte  Carlo  dose  calculation  result  and  the  dashed  line  denotes  the  conventional 
dose  calculation  result. 

FSPB  Monte  Carlo 


(a) 
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Figure  2:  These  are  the  Gaussian  distributions  for  setup  uncertainty  as  obtained  from  the  literature  for 
one  specific  treatment  site.  It  should  be  noted  that  the  distribution  is  different  along  the  x,  y  and  z 
principle  axes.  Also,  these  curves  were  generated  using  the  subroutines  from  our  Monte  Carlo 
software  so  they  serve  as  a  quality  assurance  test  of  our  code. 
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Figure  3:  This  is  a  schematic  of  our  implementation  of  setup  uncertainty  in  our  Monte  Carlo  dose 
calculation  code. 


Incident  electron 

-  initial  position 

-  initial  directions 


Electron  after  perturbation 

-  same  direction 

-  new  position  L  \ 


Model  of  setup  uncertainty 
-  Gaussian  distributions 


Figure  4:  Axial  and  sagittal  computed  tomography  slice  showing  the  Monte  Carlo  calculated  IMRT 
dose  distributions  with  (thick  line)  and  without  (thin  line)  setup  uncertainty  included.  The  target 
(metastatic  lesion  at  T1 1-T12)  is  shown  in  green  and  the  spinal  cord  in  red. 


Figure  5:  This  is  the  dose- volume  histogram  of  the  dose  distributions  shown  in  Figure  4.  There  is  an 
increase  in  dose  to  the  spinal  cord.  The  calculation  with  setup  uncertainty  is  shown  in  by  the  curves 
with  closed  circles  and  without  setup  uncertainty  by  the  solid  line. 


Dose  (Gy) 
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Figure  6:  Axial  and  sagittal  computed  tomography  slice  showing  the  Monte  Carlo  calculated  IMRT 
dose  distributions  with  (thick  line)  and  without  (thin  line)  setup  uncertainty  included.  The  target  is 
shown  in  green  and  the  spinal  cord  in  red. 


n' 


Figure  7:  This  is  the  dose- volume  histogram  of  the  dose  distributions  shown  in  Figure  6.  Note  the 
increase  in  target  dose  homogeneity.  However,  there  is  also  a  corresponding  slight  increase  in  dose  to 
the  spinal  cord.  The  calculation  with  setup  uncertainty  is  shown  in  by  the  curves  with  closed  circles 
and  without  setup  uncertainty  by  the  solid  line. 
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Abstract.  This  work  investigates  the  feasibility  of  optimizing  energy-  and  intensity-modulated 
electron  beams  for  radiation  therapy.  A  multileaf  collimator  (MLC)  specially  designed  for 
modulated  electron  radiotherapy  (MERT)  was  investigated  both  experimentally  and  by  Monte  Carlo 
simulations.  An  inverse-planning  system  based  on  Monte  Carlo  dose  calculations  was  developed  to 
optimize  electron  beam  energy  and  intensity  to  achieve  dose  conformity  for  target  volumes  near  the 
surface.  The  results  showed  that  an  MLC  with  5  mm  leaf  widths  could  produce  complex  field  shapes 
for  MERT.  Electron  intra-  and  inter-leaf  leakage  had  negligible  effects  on  the  dose  distributions 
delivered  with  the  MLC,  even  at  shallow  depths.  Focused  leaf  ends  reduced  the  electron  scattering 
contributions  to  the  dose  compared  with  straight  leaf  ends.  As  anticipated,  moving  the  MLC 
position  toward  the  patient  surface  reduced  the  penumbra  significantly.  There  were  significant 
differences  in  the  beamlet  distributions  calculated  by  an  analytic  3-D  pencil  beam  algorithm  and 
the  Monte  Carlo  method.  The  Monte  Carlo  calculated  beamlet  distributions  were  essential  to  the 
accuracy  of  the  MERT  dose  distribution  in  cases  involving  large  air  gaps,  oblique  incidence  and 
heterogeneous  treatment  targets  (at  the  tissue-bone  and  bone-lung  interfaces).  To  demonstrate  the 
potential  of  MERT  for  target  dose  coverage  and  normal  tissue  sparing  for  treatment  of  superficial 
targets,  treatment  plans  for  a  hypothetical  treatment  were  compared  using  photon  beams  and  MERT. 

(Some  figures  in  this  article  are  in  colour  only  in  the  electronic  version;  see  www.iop.org) 


1.  Intr  oduction 

Photon  beams  have  been  an  effective  modality  for  breast  cancer  treatment  in  radiation  therapy. 
Although  such  conventional  treatment  with  tangential  photon  fields  has  been  successful,  the 
following  two  problems  (or  potential  areas  of  improvement)  remain: 

(a)  The  inclusion  of  the  lung  and  other  normal  tissues,  and  sometimes  of  a  small  volume  of 
the  heart,  in  the  high-dose  volume  due  to  tumour  location,  patient  size  or  in  the  case  of 
chest-wall  treatments. 

(b)  High  exit  or  scatter  dose  to  the  normal  structures  such  as  the  lung,  the  heart  and  the 
contralateral  breast. 

Advances  in  the  state  of  the  art  of  computer-controlled  medical  linear  accelerators  have 
recently  become  available  that,  along  with  newly  developed  treatment  planning  techniques, 
may  provide  significant  improvements  in  the  delivery  and  control  of  external  beam  radiation 
through  beam-intensity  modulation  (Boesecke  et  al  1988,  Brahme  1988,  Convery  and 
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Rosenbloom  1992,  Leibel  et  al  1992,  Webb  1992,  1997,  LoSasso  et  al  1993,  Powlis  et  al 
1993,  Chui  et  al  1994,  Mageras  et  al  1994,  Brewster  et  al  1995,  Fraass  et  al  1995,  Kutcher 
et  al  1995,  Mackie  et  al  1995,  McShan  et  al  1995,  Ling  et  al  1996,  Boyer  et  al  1997).  It  is 
expected  that  using  photon  IMRT,  the  problem  (a)  above  may  be  significantly  improved  but  (b) 
may  become  more  serious  as  treatment  time  increases  with  the  number  of  fields/segments  used 
(increased  leakage  or  scattering  dose).  Using  the  modulated  electron  radiotherapy  (MERT) 
technique  (Lief  et  al  1996,  Hyodynmaa  et  al  1996,  Zackrisson  and  Karlsson  1996,  Asell  et  al 
1997,  Ebort  and  Hoban  1997,  Karlsson  et  al  1998,  1999),  on  the  other  hand,  problem  (a)  may 
also  be  significantly  improved  and  problem  (b)  may  almost  be  eliminated  due  to  the  nature  of 
the  electron  beams. 

In  the  optimization  process  of  MERT,  dose  conformity  along  the  beam  direction  can  be 
achieved  by  modulating  the  electron  incident  energy,  making  use  of  the  sharp  dose  fall-off 
feature.  A  drawback  is  its  large  penumbra  at  large  depths.  Traditionally,  electron  beams  are 
shaped  using  a  cutout  (or  blocks)  and  beam  penetration/intensity  may  be  modified  using  a  bolus. 
However,  it  is  time-consuming  to  make  such  beam  modifiers  and  the  treatment  time  would  be 
significantly  increased  if  such  beam  modifiers  are  used  for  MERT.  Efforts  have  been  made  to 
use  computer-controlled  MLC  for  electron  beam  modulation.  The  recent  results  by  Karlsson 
et  al  (1999)  showed  that  by  replacing  the  air  in  the  treatment  head  with  a  low-cost,  custom- 
made  helium  balloon,  the  beam  penumbral  width  (20/80)  was  reduced  from  18  to  1 1  mm  at 
80  cm  SSD.  The  beam  characteristics  are  affected  by  the  position  of  the  MLC.  However,  by 
replacing  the  air  between  the  MLC  and  the  patient  with  a  helium  balloon,  the  beam  penumbra 
become  almost  the  same  as  that  achieved  by  electron  beam-shaping  with  an  electron  applicator 
that  extends  to  the  patient  skin  surface  (Karlsson  et  al  1999).  This  means  that  many  of  the 
techniques  so  far  developed  with  computer-controlled  MLC  and  our  experience  with  MLC 
photon  beam  modulation  can  be  adopted  for  use  with  MERT. 

The  calculation  of  dose  distributions  for  electron  beam  radiotherapy  planning  is 
challenging  because  electron  scattering  is  strongly  affected  by  changes  in  density  and 
composition  in  the  patients.  The  3D  pencil  beam  algorithm  (Hogstrom  et  al  1981)  is  a  fast 
analytical  algorithm  which  has  been  adopted  by  most  treatment  planning  systems.  However, 
it  has  limitations  with  small  irregular  electron  fields  and  in  the  presence  of  inhomogeneities 
(Cygler  et  al  1987,  Bielajew  et  al  1987,  Mah  et  al  1989,  Mackie  et  al  1994,  Ma  et  al  1999). 
The  Monte  Carlo  simulation  has  been  demonstrated  to  be  a  viable  option  for  such  complex 
situations,  and  also  the  only  way  to  take  into  account  back-scattering  from  denser  materials  in 
a  patient  (e.g.  bone  or  metal  inserts)  (Shortt  et  al  1986,  Cygler  et  al  1987,  Mackie  et  al  1994, 
Kawrokaw  et  al  1996,  Mohan  1997,  Kapur  1999,  Ma  et  al  1999).  The  EGS4/BEAM  system 
was  developed  for  the  simulation  of  radiotherapy  beams  from  various  radiotherapy  treatment 
units,  such  medical  accelerators  (Rogers  et  al  1995).  Excellent  agreement  (1-3%)  has  been 
achieved  between  the  Monte  Carlo  dose  distributions  calculated  using  the  simulated  particle 
phase-space  data  and  measurements  (Rogers  et  al  1995,  Kapur  et  al  1998,  Zhang  et  al  1999, 
Ma  et  al  1999).  We  have  installed  a  Monte  Carlo  patient  dose  calculation  tool  on  a  clinical 
treatment  planning  system  (Ma  et  al  1999)  and  used  Monte  Carlo  for  treatment  planning  and 
dose  delivery  validation.  This  has  reduced  the  uncertainty  of  the  accelerator  output  for  small 
irregular  field  electron  beams  from  up  to  10%  to  about  3%  (Ma  et  al  1997,  Kapur  et  al  1998). 

Conformal  radiotherapy  was  initially  used  to  limit  the  normal  tissue  dose  by  conforming 
the  treatment  field  to  the  beam’s-eye-view  projection  of  the  target  volume  (Takahashi  1965). 
For  photon  beams,  the  MLC  was  used  to  collimate  the  fields  and  later  to  modulate  the  beam 
intensity  in  the  field  (Boesecke  et  al  1988,  Brahme  1988,  Convery  and  Rosenbloom  1992, 
Leibel  et  al  1992,  Webb  1992,  1997,  LoSasso  et  al  1993,  Powlis  et  al  1993,  Chui  et  al 
1994,  Mageras  et  al  1994,  Brewster  et  al  1995,  Kutcher  et  al  1995,  Mackie  et  al  1995, 
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Figure  1.  A  prototype  of  an  electron  MLC  mounted  on  the  bottom  scraper  of  a  25  cm  x  25  cm 
applicator  on  a  Varian  Clinac  2100C  accelerator.  The  MLC  has  30  pairs  of  steel  leaves  and  the  leaf 
positions  are  fixed  by  the  tightening  screws. 


McShan^/tf/  1995,  Ling  et  al  1996,  Boyer  et  al  1997).  There  have  been  a  few  studies 
on  electron  conformal  therapy.  Tailoring  dose  distributions  using  electron  beams  requires 
substantial  beam  manipulation,  due  to  their  scattering  characteristics.  Such  manipulation  is 
already  possible  with  radiation  sources  such  as  microtrons  where  preferential  energy  selection 
and  magnetically  scanned  pencil  beams  are  possible  (Lief  et  al  1996).  Both  intensity-  and 
energy/intensity-modulated  electron  beams  have  been  investigated  to  conform  the  dose  to  the 
target  near  the  surface  (Lief  et  al  1996,  Hyodynmaa  et  al  1996,  Zackrisson  and  Karlsson  1996, 
Ebert  and  Hoban  1997,  Karlsson  et  al  1998).  More  recent  work  has  studied  the  combination  of 
photon  IMRT  and  MERT  for  targets  at  greater  depth  (Karlsson  et  al  1999).  Using  the  helium- 
balloon  technique  together  with  a  computer-controlled  MLC,  it  may  be  possible  to  deliver  a 
set  of  intensity-modulated  beams  with  different  energies  and  incident  angles. 

In  this  work,  we  investigate  the  feasibility  of  optimizing  energy-  and  intensity-modulated 
electron  beams  for  radiotherapy  treatment.  We  report  here  our  Monte  Carlo  studies  of 
a  multileaf  collimator  specially  designed  for  MERT  and  some  preliminary  experimental 
results.  We  also  report  on  the  dose  calculation  algorithms  and  their  effects  on  treatment  plan 
optimization  for  MERT.  We  will  discuss  the  differences  in  the  beam  characteristics  between 
a  photon  MLC  and  an  electron  MLC.  We  will  compare  the  dose  distributions  between  a 
conventional  tangential  photon  treatment  plan  and  a  MERT  treatment  plan  for  a  hypothetical 
breast  treatment  to  demonstrate  the  potential  of  MERT  for  target  dose  coverage  and  normal 
tissue  sparing. 


2.  Materials  and  method 

2.1.  The  prototype  electron  MLC 

We  have  developed  a  prototype  manual-driven  electron  MLC  for  the  beam  delivery  for  MERT 
As  shown  in  figure  1 ,  the  electron  MLC  consists  of  30  steel  leaf  pairs,  which  were  made  from 
the  off-the-shelf  steel  bars  for  convenience  and  cost-effectiveness.  Each  leaf  is  0.476  cm  wide, 
20  cm  long  and  2.54  cm  thick  with  straight  edges  and  ends.  The  leaves  were  mounted  on  a 
steel  frame,  which  can  be  attached  to  the  bottom  scraper  of  a  25  cm  x  25  cm  electron  applicator 
on  a  Varian  Clinac  2100C  (Varian  Oncology  Systems,  Palo  Alto,  CA).  The  leaves  can  slide  in 
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the  steel  frame  and  the  leaf  positions  can  be  easily  set  using  a  pre-cut  cardboard  for  a  beam 
segment.  The  field  shape  is  maintained  by  tightening  the  screws  from  the  side.  The  largest 
radiation  field  available  using  the  electron  MLC  was  15.7  cm  x  15.7  cm  projected  at  100  cm 
source-surface  distance  (SSD). 

Because  of  the  existing  electronic  device  for  inserting  the  electron  cutout,  the  leaves  could 
not  be  placed  at  the  last  scraper  level  without  modifying  the  existing  applicator  geometry. 
Instead,  the  steel  frame  was  inserted  using  the  electron  cutout  mount  and  the  leaves  were 
placed  immediately  above  the  last  scraper.  This  resulted  in  a  slightly  greater  air  gap  (7  cm) 
between  the  bottom  of  the  electron  MLC  leaves  and  the  phantom  surface  (assuming  a  97  cm 
SSD)  compared  with  that  of  an  electron  cutout  (5  cm  for  a  100  cm  SSD).  The  projected  leaf 
width  for  a  5  cm  air  gap  is  0.5  cm,  while  for  a  7  cm  air  gap  the  projected  leaf  width  is  0.51  cm 
(e.g.  for  the  current  configuration  at  97  cm  SSD).  Further  modifications  are  needed  to  the 
electron  applicator  geometry  in  order  to  lower  the  electron  MLC  leaves.  The  ideal  location 
for  the  MLC  leaves  is  the  last  scraper  since  electron  cutouts  will  no  longer  be  needed  if  an 
electron  MLC  is  in  place. 

Film  measurement  was  performed  to  study  the  characteristics  of  the  electron  beams 
collimated  by  the  electron  MLC.  The  film  was  calibrated  following  the  AAPM  TG-25 
recommendations  (AAPM  1991)  and  the  exposures  were  taken  by  placing  film  at  different 
depths  in  a  solid  water  phantom.  The  film  was  scanned  using  a  film  scanner  which  has  a 
spatial  resolution  of  about  0.15  cm. 

2.2.  The  Monte  Carlo  beam  simulation 

We  have  used  the  EGS4  (Nelson  et  al  1985)  user  code  BEAM  for  the  accelerator  head 
simulation.  Detailed  descriptions  of  the  software  can  be  found  in  Rogers  et  al  (1995).  A 
detailed  description  of  the  clinical  implementation  of  the  Monte  Carlo  method  at  the  Stanford 
Medical  Center  was  given  in  a  previous  publication  (Ma  et  al  1999). 

For  this  work,  we  have  used  the  previously  simulated  Monte  Carlo  beam  data  for  6,  12 
and  20  MeV  electron  beams  from  a  Varian  Clinac  2100C  linear  accelerator  and  for  6  MV 
photon  beams  from  a  2300CD  accelerator  (Varian  Oncology  Systems,  Palo  Alto,  CA).  The 
dimensions  and  materials  for  the  accelerator  components  were  incorporated  according  to  the 
manufacturer’s  specifications.  Electron  beams  emerging  from  the  vacuum  exit  window  were 
assumed  to  be  monoenergetic  and  monodirectional  with  a  beam  radius  of  0. 1-0.2  cm  (Kapur 
et  al  1998).  The  energy  cutoffs  for  electron  transport  in  the  accelerator  simulation  (ECUT 
and  AE)  were  700  keV  (kinetic  +  rest  mass)  and  for  photon  transport  (PCUT  and  AP)  10  keV. 
The  electron  transport  step  length  was  confined  such  that  the  maximum  fractional  energy  loss 
per  electron  step  is  4%  (i.e.  ESTEPE  =  0.04).  The  ICRU  recommended  compositions  and 
stopping  power  values  were  used  for  the  materials  in  the  accelerator  simulations  (ICRU  1984). 
The  phase-space  data  were  scored  at  a  plane  either  immediately  above  the  photon  MLC  or 
above  the  lowest  scraper.  The  number  of  particles  was  about  2-30  million  in  an  electron  beam 
file  and  about  50  million  in  a  photon  file. 

Field  shaping  by  the  photon  MLC  or  electron  MLC  was  further  simulated  using  the  BEAM 
component  module  MLC.  MLC  could  simulate  either  straight  or  ‘double  focused’  leaf  edges 
and  ends.  In  this  work,  we  have  simulated  electron  beams  collimated  by  a  photon  MLC  with 
both  straight  and  double  focused  MLC  leaf  shapes.  The  leaves  were  7.5  cm  thick  and  made  of 
tungsten.  The  leaf  center  was  49  cm  from  the  isocentre.  The  intervening  air  in  the  accelerator 
and  between  the  MLC  and  the  isocentre  was  in  some  cases  replaced  with  helium  to  investigate 
the  effect  of  electron  scattering  in  the  air.  In  the  simulations  of  the  electron  beams  collimated  by 
an  electron  MLC,  the  leaves  were  placed  on  the  bottom  scraper  of  a  25  cm  x  25  cm  applicator 
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with  a  7  cm  air  gap  between  the  bottom  surface  of  the  MLC  and  the  isocentre.  Tungsten  leaves 
of  1 .5  cm  thickness  with  straight  edges  and  ends  were  used  in  ail  the  simulations  and  the  phase 
space  data  were  used  in  the  subsequent  dose  calculations  except  for  the  leaf  leakage  study 
where  different  leaf  materials  and  thicknesses  were  investigated  for  the  electron  MLC. 


2.5.  The  Monte  Carlo  dose  calculation 

The  EGS4  user  code,  MCDOSE  (Ma  etal  1 999),  was  used  in  this  work  for  the  dose  calculations. 
MCDOSE  was  designed  for  dose  calculations  in  a  3D  rectilinear  voxel  geometry.  Voxel 
dimensions  were  completely  variable  in  all  three  directions.  Every  voxel  (volume  element) 
could  be  assigned  to  a  different  material.  The  cross-section  data  for  the  materials  used  were 
available  in  a  pre-processed  PEGS4  cross-section  data  file.  The  mass  density  of  the  material 
in  a  MCDOSE  calculation  was  varied  based  on  the  patient's  CT  data  although  the  density 
effect  corrections  for  the  stopping  powers  of  the  material  remain  unchanged  (Ma  et  al  1999). 
The  voxel  dimensions  and  materials  were  defined  in  a  MCDOSE  input  file  together  with  the 
transport  parameters  such  as  the  energy  cutoffs  (ECUT  and  PCUT),  the  maximum  fractional 
energy  loss  per  electron  step  (ESTEPE),  and  the  parameters  required  by  PRESTA  (Bielajew 
and  Rogers  1987).  Several  variance  reduction  techniques  have  been  implemented  in  the 
MCDOSE  code  to  improve  the  calculation  efficiency.  These  include  photon  interaction  forcing, 
particle  splitting,  Russian  roulette,  electron  range  rejection  and  region  rejection,  particle  track 
displacement  and  rotation,  and  correlated  sampling.  Detailed  descriptions  of  these  techniques 
have  been  given  elsewhere  (Rogers  and  Bielajew  1990,  Ma  and  Nahum  1993,  Rogers  et  al 
1995,  Kawrakow  et  al  1996,  Keall  and  Hoban  1996,  Ma  et  al  1999). 

For  patient  dose  calculations,  the  simulation  phantom  was  built  from  the  patient’s  CT  data 
with  up  to  128  x  128  x  128  voxels  (uniform  in  any  dimension).  The  side  of  a  voxel  varied 
from  0.2  to  0.4  cm.  A  separate  program  was  developed  to  convert  the  patient’s  CT  data  from 
the  FOCUS  treatment  planning  system  (Computerized  Medical  Systems,  St  Louis,  MO)  to 
desired  dimensions,  material  types  and  densities.  The  organ  contours  were  also  obtained  for 
dose  calculation  and  analysis.  The  phase-space  data  obtained  from  a  BEAM  simulation  were 
used  as  a  source  input  with  variable  source  positions  and  beam  incident  angles.  To  simulate 
the  dose  distribution  of  a  finite  size  beamlet  used  by  the  inverse  planning  process,  particles 
were  transported  to  the  MLC  plane  and  only  those  within  the  beamlet  area  (=  1  cm  x  1  cm 
projected  at  100  cm  SSD)  were  allowed  to  go  through.  This  ignored  the  bremsstrahlung 
photon  leakage  and  electron  scattering  by  the  leaf  ends  in  the  optimization  process  (the  effect 
was  corrected  in  the  final  dose  calculation,  as  discussed  below).  After  optimization,  a  leaf 
sequence  was  generated  using  a  modified  ‘step  and  shoot’  algorithm  based  on  our  early  work 
(Ma  et  al  1998).  The  final  MERT  dose  distribution  was  computed  based  on  an  intensity  map 
(a  2D  distribution  of  particle  weighting  factors)  reconstructed  from  the  leaf  sequence.  The 
bremsstrahlung  leaf  leakage  effect  was  included  in  the  intensity  map  using  the  leaf  sequence 
and  pre-calculated  leaf  leakage  data  for  1.5  cm  thick  tungsten  leaves.  MCDOSE  produced 
data  files  that  contained  geometry  specifications  such  as  the  number  of  voxels  in  all  the  three 
directions  and  their  boundaries  as  well  as  the  dose  values  and  the  associated  (lcr)  statistical 
uncertainties  in  the  individual  voxels  and  organs  (structures).  The  EGS4  transport  parameters 
were  ECUT  =  AE  =  700  keV,  PCUT  =  AP  =  10  keV  and  ESTEPE  =  0.04.  The  number  of 
particle  histories  simulated  ranged  from  2  million  to  30  million  for  a  MERT  treatment.  The 
la  statistical  uncertainty  in  the  dose  was  generally  2%  or  smaller  of  the  £>max  value.  The  CPU 
time  required  for  a  MERT  simulation  was  about  1-3  h  on  a  Pentium  III  450  MHz  PC  with  the 
variance  reduction  option  switched  on. 
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A  B 

Figured  Beam  intensity  distributions  measured  by  film  on  the  surface  of  a  solid  water  phantom 
for  6,  12  and  20  MeV  electrons.  The  MLC  leaf  positions  for  the  electron  fields  are  also  shown 
(bottom  left). 


2.4.  The  optimization  process 

The  treatment  planning  optimization  system  used  in  this  work  is  a  home-developed  system 
based  on  the  work  by  Jiang  (1998).  First,  the  planner  inputs  the  patient  geometry  and  defines 
the  treatment  setup,  such  as  the  beam  energy,  number  and  orientations  of  beams,  etc.  The 
target  volume  and  the  critical  structures  are  defined  by  the  clinician.  A  reference  monitor  unit  is 
assigned  to  each  open  rectangular  beam  and  the  dose  deposition  coefficients,  which  are  defined 
as  the  dose  contribution  from  a  beamlet  to  a  point,  are  calculated  using  the  MCDOSE  code. 

Second,  using  the  calculated  dose  deposition  coefficients  as  input,  the  optimal  intensity 
profile  for  each  beam  is  achieved  using  a  gradient  method  to  minimize  the  objective  function. 
For  the  target  area,  a  quadratic  form  of  objective  function  is  specified.  In  addition,  two 
target  dose-uniformity  constraints  are  used  to  ensure  a  uniform  target  dose  distribution  and 
to  distinguish  the  clinical  importance  of  cold  and  hot  spots.  For  the  critical  structures, 
maximum-dose  constraint  and  several  levels  of  dose- volume  constraints  are  assigned  to  each 
structure.  For  each  objective  function  and  constraint,  an  importance  weight  relative  to  the 
target  objective  function  is  assigned.  All  the  constraints  are  mathematically  transformed  to 
the  penalty  functions  of  quadratic  forms.  The  augmented  objective  function,  which  should  be 
minimized,  is  a  combination  of  the  original  objective  functions  and  all  penalty  functions.  The 
results  of  the  optimization  process  are  the  intensity  profiles  for  the  individual  fields  (different 
incident  energies  and  gantry  angles).  The  same  optimizer  has  been  used  for  photon  beam 
optimization  with  the  Monte  Carlo  method  and  a  finite-size  pencil  beam  algorithm  (Jiang 
1998,  Jiang  et  al  1999,  Pawlicki  et  al  1999). 

3.  Resultsand  discussion 

3.1.  Characteristics  of  electron  beams  collimated  by  an  electron  MLC 

Figure  2  shows  the  electron  fields  collimated  by  the  prototype  electron  MLC  for  6,  12  and 
20  MeV  electron  beams  on  a  Varian  Clinac  2 1 00C  machine.  For  convenience,  a  photo  showing 
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Figure  3.  Beam  profiles  measured  by  film  on  the  surface  of  a  solid  water  phantom  for  the  6  MeV 
and  20  MeV  electron  fields  shown  in  figure  2:  (a)  20  MeV  along  A-A;  ( b )  20  MeV  along  B-B; 
(c)  6  MeV  along  A-A;  (d)  6  MeV  along  B-B. 


the  MLC  leaf  positions  for  the  field  shape  is  also  included  in  figure  2.  Figure  2  shows  the 
film  measurement  at  the  surface  of  a  solid  water  phantom  (97  cm  SSD)  for  6,  12  and  20  MeV 
electron  beams.  Figure  3  shows  the  measured  profiles  on  the  phantom  surface  along  A-A  and 
B-B  for  the  6  MeV  and  20  MeV  electron  fields  shown  in  figure  2.  Figure  4  shows  the  beam 
profiles  at  2  cm  depth  in  the  solid  water  phantom.  It  can  be  seen  that  for  a  20  MeV  electron 
beam,  0.5  cm  leaf  shapes  are  still  distinguishable  on  the  surface  but  become  very  blurred  at 
2  cm  depth.  For  a  6  MeV  electron  beam,  however,  the  effect  of  electron  scattering  becomes  so 
severe  that  a  leaf  width  smaller  than  1.0  cm  will  not  result  in  any  improvement  in  the  spatial 
resolution.  However,  a  small  leaf  width  may  have  the  advantage  of  defining  the  field  more 
precisely  in  the  direction  perpendicular  to  the  leaves. 

Based  on  these  experimental  results,  we  further  performed  Monte  Carlo  simulations 
of  electron  fields  collimated  by  1  cm  wide  leaves  to  study  the  effect  of  material  type  and 
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Figure4.  Beam  profiles  measured  at  2  cm  depth  in  a  solid  water  phantom  for  the  6MeV  and  20  MeV 
electron  fields  shown  in  figure  2:  {a)  20  MeV  along  A-A;  ( b )  20  MeV  along  B-B;  (c)  6  MeV  along 
A-A;  and  (d)  6  MeV  along  B-B. 


leaf  thickness.  Although  the  beam  penumbral  widths  did  not  change  significantly  for  leaf 
thicknesses  smaller  than  2  cm,  the  beam  intensity  outside  the  field  was  affected  by  the  leaf 
thickness  and  the  atomic  number  of  the  leaf  material.  As  shown  in  figure  5  for  a  20  MeV  electron 
beam,  1 .5  cm  thick  zinc  reduced  the  electron  fluence  outside  the  field  to  about  5%  of  the  central 
axis  value  (figure  5(a)).  These  electrons  were  mainly  generated  by  the  bremsstrahlung  photons 
in  the  MLC  leaves.  This  was  confirmed  by  the  photon  fluence  as  shown  in  figure  5(&),  where 
1 .5  cm  zinc  MLC  leaves  resulted  in  about  60%  higher  photon  fluence  outside  the  field  compared 
with  the  central  axis  photon  fluence.  Some  electrons  were  also  scattered  off  the  leaf  ends  and 
by  air.  For  1.5  cm  copper,  1.5  cm  lead  and  2  cm  steel,  the  electron  fluence  was  about  2.5%  of 
the  central  axis  value.  The  electron  fluence  was  reduced  to  about  1 .5%  if  the  leaves  were  made 
of  1.5  cm  tungsten.  This  was  reflected  by  the  30%  smaller  photon  fluence  under  the  tungsten 
MLC  leaves  compared  with  the  central  axis  photon  fluence.  Clearly,  tungsten  was  superior  to 


Modulated  electron  beam  for  radiotherapy 


2301 


Clinac  2100C  20  MeV  10cm  x  10cm  field 


(a) 


Clinac  2100C  20  MeV  10cm  x  10cm  field 


(*) 

Figure  5.  Monte  Carlo  simulated  electron  ( a )  and  photon  ( b )  planar  fluence  in  the  penumbral 
region  and  outside  the  treatment  field  fora  Varian  Clinac  2100C  20 MeV  electron  beam  collimated 
by  an  electron  MLC  of  different  leaf  materials  and  thicknesses.  The  air  gap  between  the  electron 
MLC  and  the  scoring  plane  is  7  cm. 


other  materials  in  terms  of  leaf  leakage.  If  we  increased  the  tungsten  leaf  thickness  to  2  cm 
the  electron  fluence  would  be  reduced  to  less  than  1%  of  the  central  axis  value  and  the  photon 
leakage  would  be  reduced  to  about  50%  of  the  central  axis  value  (not  shown). 

To  study  the  overall  effect  of  the  leaf  leakage,  leaf  scattering,  air  scattering  and  the  extended 
source  in  an  electron  beam,  we  compared  the  dose  distributions  for  single  fields  and  multiple 
abutting  fields  collimated  by  an  electron  MLC  with  1.5  cm  thick  tungsten  leaves.  Figure  6 
shows  the  Monte  Carlo  calculated  dose  distributions  for  a  single  4  cm  x  4  cm  electron  field 
and  a  multiple  abutting  field  of  the  same  size  formed  by  four  1  cm  x  4  cm  electron  fields.  For 
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Clinac  2100C  20MeV,  7cm  gap 


Clinac  2100C  20  MeV,  7cm  gap 


Clinac  2100C  6MeV,  7cm  gap 


(C) 


Figure  6.  Monte  Carlo  calculated  dose  distributions  in  a  water  phantom  for  Varian  Clinac  6  and 
20  MeV  electron  beams  collimated  by  an  electron  MLC  of  1 .5  cm  thick  tungsten  leaves  for  a  single 
4  cm  x  4  cm  electron  field  and  a  4  cm  x  4  cm  field  formed  by  four  1  cm  x  4  cm  electron  fields: 
(a)  dose  at  surface  for  a  20  MeV  beam;  ( b )  dose  at  3  cm  depth  for  a  20  MeV  beam;  (c)  dose  at 
surface  for  a  6  MeV  beam. 


a  20  MeV  electron  beam,  the  dose  at  the  phantom  surface  for  the  abutting  field  shows  about 
4%  fluctuation  compared  with  a  single  electron  field  (figure  6(a)).  This  is  potentially  due  to 
the  effect  of  leaf  shape  and  extended  source.  The  dose  outside  the  field  for  the  abutting  field  is 
about  three  times  higher  than  that  for  the  single  field,  which  is  mainly  caused  by  the  leaf  leakage 
due  to  the  longer  beam-on  time  to  deliver  the  four  1  cm  x  4  cm  fields  and  electron  scattering 
off  the  leaf  ends.  The  dose  at  3  cm  depth  shows  little  difference  between  the  abutting  field  and 
the  single  field  except  for  the  dose  near  the  field  edges  and  outside  the  field  (figure  6(b)).  For  a 
6  MeV  electron  beam,  the  dose  at  the  phantom  surface  for  the  abutting  field  is  almost  the  same 
as  that  for  the  single  field  (figure  6(c)).  The  dose  outside  the  field  for  the  abutting  field  is  only 
slightly  higher  than  that  for  the  single  field.  The  effect  of  leaf  leakage  is  very  small  for  a  6  MeV 
beam  and  the  dose  immediately  outside  the  field  is  thought  to  be  mainly  due  to  the  effect  of 
electron  scattering  in  the  air.  It  seems  that  field  abutting  with  1  cm  beamlets  collimated  by 
an  electron  MLC  can  provide  adequate  beam  characteristics  for  MERT  for  the  beam  energies 
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investigated.  However,  the  dose  outside  the  field  needs  to  be  minimized  through  beam  energy 
and  leaf  sequence  optimization. 

3.2.  Comparisons  of  a  photon  MLC  and  an  electron  MLC 

There  have  been  studies  on  electron  beam  collimation  using  a  photon  MLC  (Karlsson  et  al 
1999).  One  of  the  advantages  of  using  a  photon  MLC  is  the  possibility  of  combining  both 
photon  and  electron  beams  in  the  same  plan.  An  essential  requirement  for  matching  a  photon 
beam  and  an  electron  beam  at  different  depths  is  that  both  beams  share  the  same  source 
position.  Karlsson  et  al  (1999)  proposed  several  modifications  to  the  design  of  a  Varian 
Clinac  2300CD  accelerator,  one  of  which  was  to  replace  the  intervening  air  with  helium.  This 
could  significantly  reduce  the  effect  of  electron  scattering  in  the  air  on  the  beam  penumbra. 
However,  filling  the  accelerator  head  with  helium  requires  major  modifications  to  the  existing 
accelerator  design.  In  this  work,  we  have  investigated  an  alternative  solution — a  thin  leaf 
MLC  at  the  electron  cutout  level  to  reduce  the  air  scattering  effect.  As  can  be  seen  in  figure  7, 
the  unfocused  MLC  leaf  ends  could  scatter  the  electrons  very  significantly  to  degrade  the 
beam  characteristics  near  the  field  edges.  The  Varian  MLC  has  rounded  leaf  ends,  which 
are  expected  to  have  similar  dosimetric  characteristics  as  the  unfocused  MLC  studied  here. 
Focused  leaf  ends  could  greatly  improve  the  beam  edges  and  provided  even  slightly  better 
dose  profiles  inside  the  field  for  a  20  MeV  electron  beam  compared  with  an  electron  MLC 
(figures  7(a)-(c)),  primarily  due  to  the  reduction  of  electron  scattering  in  the  accelerator  head 
(helium  versus  air).  The  dose  outside  the  field  was  slightly  lower  for  the  electron  MLC  than  for 
the  photon  MLC.  For  a  6  MeV  beam,  an  electron  MLC  gave  slightly  better  surface  dose  profiles 
both  inside  and  outside  the  field  than  the  focused  and  unfocused  photon  MLC.  However,  the 
dose  profiles  became  practically  similar  at  the  depth  of  the  maximum  dose  and  greater  depths 
(not  shown).  Note  that  in  these  comparisons,  we  have  placed  the  phantom  surface  at  20  cm 
below  the  photon  MLC  and  7  cm  below  the  electron  MLC  to  minimize  the  effect  of  electron 
scattering  in  the  air  or  helium  between  the  MLC  and  the  phantom.  It  is  evident  that  an  electron 
MLC  will  have  similar  dosimetric  characteristics  as  a  photon  MLC  with  focused  leaf  ends  but 
without  the  need  to  replace  the  air  in  the  accelerator  head  with  helium. 

3.3.  Comparison  of  beamlet  distributions 

The  accuracy  of  the  beamlet  distribution  calculation  may  play  an  important  role  in  the  treatment 
planning  optimization  process.  Ma  et  al  (1999)  reported  significant  differences  in  the  final  dose 
distributions  of  the  optimized  treatment  plans  computed  by  a  commercial  inverse  treatment 
planning  system  with  a  finite-size  pencil  beam  and  the  Monte  Carlo  method.  Pawlicki  et  al 
(1999)  demonstrated  that  inaccurate  beamlet  distributions  may  result  in  under-dosing  in  the 
target  and  over-dosing  in  the  adjacent  critical  structures,  and  using  the  Monte  Carlo  calculated 
beamlets  could  potentially  reduce  the  uncertainty  in  the  photon  IMRT  dose  distributions.  This 
was  demonstrated  again  by  Jeraj  and  Keall  (1999)  using  a  Monte  Carlo  dose  calculation  based 
inverse  planning  algorithm. 

It  has  been  shown  that  the  electron  beam  dose  distributions  calculated  by  the  pencil  beam 
algorithm  as  implemented  in  some  commercial  treatment  planning  systems  could  be  fairly 
uncertain  in  the  regions  near  material  interfaces  and  inhomogeneities  (Cygler  et  al  1987, 
Shortt  et  al  1986,  Mackie  et  al  1994,  Ma  et  al  1999).  We  have  computed  the  beamlet  dose 
distributions  using  the  3D  pencil  beam  as  implemented  in  the  FOCUS  treatment  planning 
system  (Computerized  Medical  Systems,  St  Louis,  MO)  and  compared  them  with  the  Monte 
Carlo  calculated  beamlets.  Figure  8  shows  the  dose  distributions  calculated  using  the  Monte 
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Figure7.  Monte  Carlo  calculated  dose  distributions  in  a  water  phantom  for  a  10cm  x  10cm  field 
collimated  by  an  electron  MLC  with  1.5  cm  thick  tungsten  leaves  and  a  photon  MLC  with  7.5  cm 
thick  leaves  on  a  Varian  Clinac  2100C  accelerator:  (a)  surface  dose  for  a  20  MeV  beam;  ( b )  dose 
at  3  cm  depth  for  a  20  MeV  beam;  (c)  dose  at  6  cm  depth  for  a  20  MeV  beam;  ( d )  surface  dose  for 
a  6  MeV  beam.  The  electron  MLC  has  straight  leaf  ends.  The  photon  MLC  has  either  straight  or 
double-focused  leaves. 


Carlo  method  (a,  c,  e)  and  the  FOCUS  3D  pencil  beam  algorithm  (bt  dy  f)  for  a  1  cm  x  1  cm 
12  MeV  beamlet  incident  on  a  patient  phantom  built  from  CT  data.  For  beamlets  with  normal 
incidence  (figures  8(a)  and  {b))y  the  difference  in  the  dose  distributions  in  the  heart  was  evident: 
the  Monte  Carlo  calculated  isodose  lines  varied  with  the  heart  contour  while  the  pencil  beam 
isodose  lines  remained  symmetrical  despite  the  change  in  material  densities.  Figures  8(c) 
and  ( d )  show  the  beamlet  distributions  with  a  10  cm  air  gap.  The  difference  is  clearly  seen 
near  the  surface.  The  beamlet  distributions  again  differed  significantly  in  the  lung  for  oblique 
incidence  (figures  8(c)  and  (f )).  The  axis  of  the  beamlet  was  intentionally  placed  to  go  through 
soft  tissues  and  bones.  The  pencil  beam  isodose  lines  seemed  to  stretch  according  to  the  beam 
axis  pathlength  and  showed  no  signs  of  electron  build-down  near  the  low-density  material. 
These  results  provided  enough  evidence  to  show  that  to  ensure  the  accuracy  of  the  optimized 
dose  distributions  for  MERT  we  should  use  the  Monte  Carlo  method  to  compute  the  electron 
beamlets  for  the  inverse  planning  process. 
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Figure  8.  Dose  distributions  calculated  using  the  Monte  Carlo  method  ( a ,  c,  e)  and  the  FOCUS 
3D  pencil  beam  algorithm  (b,  d ,  /)  for  a  1  cm  x  1  cm  12  MeV  beamlet  with  normal  incidence 
{a  and  b ),  normal  incidence  plus  10  cm  air  gap  (c  and  d),  and  oblique  incidence  ( e  and  /).  The 
beamlet  size  is  defined  at  100  cm  SSD.  The  isodose  lines  shown  are  10,  20,  30,  40,  50,  60,  70,  80 
and  90%  of  the  maximum  dose  respectively. 


3.4.  MERT  versus  photons:  a  hypothetical  treatment  plan 

Modulated  electron  radiotherapy  is  a  general  purpose  technique  that  should  be  advantageous 
in  many  clinical  situations.  An  exhaustive  investigation  of  the  specific  advantages  of  MERT 
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Figure  .  Treatment  plans  for  a  hypothetical  breast  case  using  tangential  6  MV  photon  beams  ( a ) 
and  MERT  with  6, 12  and  20  MeV  electron  beams  ( b ).  Both  plans  were  calculated  using  the  Monte 
Carlo  method.  The  isodose  lines  (90,  70,  50  and  30%)  represented  50,  38.9,  25  and  16.7  Gy. 


over  traditional  treatment  modalities  on  a  site  by  site  basis  is  outside  the  scope  of  this  work. 
However,  to  demonstrate  the  possibility  of  improving  the  dose  homogeneity  in  the  target  and  the 
reduction  of  the  dose  in  surrounding  normal  tissues,  we  compare  the  dose  distributions  to  treat 
a  hypothetical  target  using  tangential  photon  beams  and  MERT.  The  purpose  of  the  comparison 
was  to  illustrate  the  concept  of  MERT  but  not  to  draw  specific  conclusions  on  the  use  of  either 
technique.  Previous  investigators  have  used  artificial  phantoms  and  hypothetical  targets  to 
mimic  different  treatment  sites  (e.g.  Hyodynmaa  etal  1996,  Asell  et  al  1997,  Ebert  and  Hoban 
1997).  We  considered  it  to  be  clinically  relevant  to  use  a  more  realistic  patient  geometry  (built 
from  CT  data)  in  our  comparison,  although  the  target  definition  and  beam  setup  are  somewhat 
arbitrary.  Figure  9  shows  the  hypothetical  treatment  plan  using  tangential  6  MV  photon  beams 
and  MERT  with  normally  incident  6,  12  and  20  MeV  electron  beams.  The  intensity  maps  for 
each  electron  beam  energy  are  shown  in  figure  10.  The  beamlet  size  was  1  cm  x  1  cm  at  100  cm 
at  isocentre.  It  is  worth  noting  here,  that  as  a  matter  of  practicality,  it  is  impossible  to  create 
the  intensity  maps  shown  in  figure  10  using  the  conventional  electron  cutout  approach  but  the 
electron  MLC  is  a  viable  alternative.  The  dose  distributions  for  both  plans  were  calculated 
using  the  Monte  Carlo  method.  The  isodose  lines  (90,  70,  50  and  30%)  were  normalized  in 
such  a  way  that  the  90%  isodose  surface  would  receive  the  prescribed  target  dose  of  50  Gy. 
For  the  photon  plan,  the  90%  dose  line  also  included  a  margin  in  the  lung  to  account  for  the 
effect  of  patient  breathing.  This  was  not  needed  for  the  electron  plan  as  the  electron  beams 
were  incident  enface  and  the  electron  beamlet  dose  distributions  do  not  vary  significantly  with 
breathing.  Figure  1 1  shows  the  dose  volume  histograms  (DVH)  for  the  hypothetical  treatment 
plans  shown  in  figure  9.  The  target  DVH  together  with  the  right  lung  DVH  are  shown  in 
figure  11(a)  (as  percentage  volume)  and  the  right  lung  DHV  and  the  ‘total  body’  (including 
everything  inside  the  external  contour)  DVH  are  shown  in  figure  1 1(6)  (as  absolute  volume). 

It  is  clear  that  MERT  provided  better  dose  homogeneity  in  the  target  region  than  tangential 
photon  beams.  Tangential  photon  beams  produced  hot  spots  in  the  target  and  cold  spots  near 
the  skin  (figure  11(a)).  MERT  significantly  reduced  the  dose  to  the  lung  relative  to  tangential 
photon  beams;  the  maximum  dose  to  the  lung  was  reduced  from  50  Gy  for  a  tangential  treatment 
to  35  Gy  for  MERT  (figure  11(a)).  However,  MERT  increased  the  volume  of  the  lung  that 
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Figure  10.  Intensity  maps  for  the  three  electron  beam  energies.  The  beamlet  size  for  each  port 
was  1  x  1  cm2.  Darker  beamlets  indicate  a  higher  weight  than  the  lighter  beamlets  and  the  grey 
scale  for  all  three  maps  is  in  absolute  terms. 


received  a  lower  dose  (10%  more  volume  received  5  Gy  and  20%  more  volume  received  2  Gy) 
compared  with  tangential  photon  beams.  The  clinical  significance  of  the  increased  lung  volume 
receiving  such  a  low  dose  needs  to  be  investigated.  On  the  other  hand,  over  150  cm3  of  lung 
received  much  less  dose  with  MERT  compared  with  tangential  photon  beams,  which  could 
result  in  reduced  lung  complications  (figure  11  (fc)).  Another  clear  benefit  with  MERT  is  the 
exclusion  of  the  surrounding  normal  tissue  from  the  high  dose  volume  (figure  11(&)).  Over 
1000  cm3  of  normal  tissue  received  10-30  Gy  less  dose  in  a  MERT  plan  compared  with  this 
tangential  photon  beam  plan. 

4.  Conclusions 

In  this  work,  we  have  investigated  the  feasibility  of  modulating  both  energy  and  intensity  of 
electron  beams  for  radiotherapy.  This  was  achieved  by  combining  electron  beams  of  different 
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(b) 

Figure  1 1 .  Dose  volume  histograms  (DVH)  for  the  breast  treatment  plans  shown  in  figure  9: 
(a)  DVH  shown  as  percentage  volume  for  the  target  (PTV)  and  the  right  lung  and  ( b )  DVH  shown 
as  absolute  volume  for  the  right  lung  and  the  ‘total  body’  which  includes  everything  inside  the 
external  contour 


nominal  energies  and  variable  intensity  distributions.  A  prototype  electron  MLC  was  built  to 
study  the  characteristics  of  MLC-collimated  electron  beams  and  the  Monte  Carlo  simulations 
were  used  to  investigate  the  effect  of  MLC  leaf  material,  thickness,  shape  and  location.  The 
beamlet  distributions  calculated  using  a  3D  electron  pencil  beam  algorithm  as  implemented 
in  a  commercial  treatment  planning  system  and  the  Monte  Carlo  method  were  compared  for 
electron  beams  of  different  energies,  extended  air  gaps,  oblique  incidence  and  heterogeneous 
geometries.  A  hypothetical  breast  case  was  used  to  compare  the  dose  distributions  using 
tangential  photons  and  MERT  for  target  coverage  (dose  homogeneity)  and  normal  tissue  sparing 
(dose  reduction  in  the  lung  and  other  surrounding  normal  tissues). 
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Our  results  showed  that  an  electron  MLC  at  the  electron  cutout  location  can  provide 
adequate  beam  collimation  for  MERT  without  the  need  to  replace  the  air  in  the  accelerator 
head  and  between  the  MLC  and  the  phantom  with  helium.  The  beam  characteristics  collimated 
by  an  electron  MLC  are  comparable  to  those  collimated  by  a  focused  photon  MLC.  However, 
the  latter  requires  the  accelerator  head  and  between  the  MLC  and  the  phantom  to  be  filled  with 
helium,  which  may  be  impractical  for  some  accelerators  because  of  the  major  modifications 
needed  to  the  structure  design.  An  electron  MLC  can  also  be  used  in  place  of  a  cutout.  The 
Monte  Carlo  method  can  accurately  simulate  particle  transport  in  cases  involving  extended  air 
gaps,  oblique  incidence  and  heterogeneous  anatomy,  and  is  therefore  suitable  for  the  beamlet 
calculation  for  MERT  treatment  optimization.  Our  preliminary  results  based  on  a  hypothetical 
breast  case  demonstrated  the  potential  of  MERT  for  uniform  target  coverage  and  normal  tissue 
sparing.  To  fully  explore  the  potential  of  MERT,  further  studies  need  to  be  carried  out  for 
realistic  clinical  cases  and  for  other  treatment  sites  such  as  the  head  and  neck. 
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Abstract.  A  new  EGS4/PRESTA  Monte  Carlo  user  code,  MCDOSE,  has  been  developed  as 
a  routine  dose  calculation  tool  for  radiotherapy  treatment  planning.  It  is  suitable  for  both 
conventional  and  intensity  modulated  radiation  therapy.  Two  important  features  of  MCDOSE 
are  the  inclusion  of  beam  modifiers  in  the  patient  simulation  and  the  implementation  of  several 
variance  reduction  techniques.  Before  this  tool  can  be  used  reliably  for  clinical  dose  calculation,  it 
must  be  properly  validated.  The  validation  for  beam  modifiers  has  been  performed  by  comparing 
the  dose  distributions  calculated  by  MCDOSE  and  the  well-benchmarked  EGS4  user  codes  BEAM 
and  DOSXYZ.  Various  beam  modifiers  were  simulated.  Good  agreement  in  the  dose  distributions 
was  observed.  The  differences  in  electron  cutout  factors  between  the  results  of  MCDOSE  and 
measurements  were  within  2%.  The  accuracy  of  MCDOSE  with  various  variance  reduction 
techniques  was  tested  by  comparing  the  dose  distributions  in  different  inhomogeneous  phantoms 
with  those  calculated  by  DOSXYZ  without  variance  reduction.  The  agreement  was  within  1.0%. 
Our  results  demonstrate  that  MCDOSE  is  accurate  and  efficient  for  routine  dose  calculation  in 
radiotherapy  treatment  planning,  with  or  without  beam  modifiers. 

(Some  figures  in  this  article  are  in  colour  only  in  the  electronic  version;  see  www .  iop .  org) 


1.  Introduction 

With  the  rapid  development  of  computer  technology,  the  Monte  Carlo  technique,  currently 
the  most  accurate  method  for  dose  calculation  (Shortt  et  al  1986,  Mackie  1990,  Rogers  and 
Bielajew  1990,  Andreo  1991,  Mackie  et  al  1994,  DeMarco  et  al  1998,  Ma  et  al  1999a),  is 
becoming  more  practical  for  use  in  radiation  therapy  treatment  planning  systems.  It  has  been 
argued  that  the  outcome  of  radiation  therapy  treatment  may  be  improved  by  using  Monte  Carlo 
dose  calculation  (Mohan  1997,  Nahum  1997),  especially  for  intensity  modulated  radiotherapy 
(IMRT)  treatment  planning  (Jeraj  and  Keall  1999,  Ma  et  al  1999a,  2000a,  b).  The  Monte  Carlo 
method  is  being  developed  for  treatment  planning  dose  calculations  by  several  groups  (Mackie 
1990,  Hartmann-Siantar  et  al  1997,  Ma  et  al  1997,  1999a,  DeMarco  et  al  1998,  Wang  et  al 
1998,  Mubata  et  al  1998,  Libby  et  al  1998,  Faddegon  et  al  1998). 

At  present,  there  are  several  general-purpose  Monte  Carlo  codes  in  widespread  use 
for  radiation  transport  simulation,  e.g.  Electron  Gamma  Shower  version  4  (EGS4)  (Nelson 
et  al  1985),  ETR AN/ITS  (Berger  and  Seltzer  1987)  and  Monte  Carlo  Af-particle  (MCNP) 
(Hendricks  1994)  and  PENELOPE  (Salvat  et  al  1996).  It  is  well  known  that  the  EGS4  code 
system  is  very  well  documented  and  it  has  been  thoroughly  benchmarked  in  the  energy  region 
of  dosimetric  interest  (Rogers  1984a,  b,  Rogers  and  Bielajew  1984,  1986,  1990).  DOSXYZ  is 
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an  EGS4  based  Monte  Carlo  simulation  code  for  calculation  dose  distribution  in  a  rectilinear 
voxel  phantom  (Rogers  et  al  1995,  Ma  et  al  1995).  It  is  also  well  benchmarked.  Comparisons 
of  the  DOSXYZ  results  with  measurements  have  been  reported  previously  (Rogers  et  al  1995, 
Kapur  et  al  1998,  Ma  1998,  Zhang  et  al  1999,  Ma  et  al  1999a,  2000a).  Unfortunately,  these 
codes  are  too  slow  to  be  acceptable  for  routine  dose  calculation  in  treatment  planning  systems 
with  existing  hardware.  Some  codes  based  on  the  Monte  Carlo  method  have  been  developed 
to  speed  up  the  dose  calculation,  such  as  Macro  Monte  Carlo  (MMC)  (Neuenschwander  and 
Born  1992,  Neuenschwander  etal  1995),  Voxel  Monte  Carlo  (VMC  and  XVMC)  (Kawrakow 
et  al  1996,  Fippel  1999)  and  Super  Monte  Carlo  (SMC)  (Keall  and  Hoban  1996a).  Further 
verification  and  improvements  are  needed  for  routine  treatment  planning  dose  calculations. 

A  new  EGS4/PRESTA  (Bielajew  and  Rogers  1987)  user  code,  MCDOSE,  has  been 
developed  for  routine  clinical  dose  calculation  (Ma  et  al  1999b,  2000c).  MCDOSE  was 
designed  as  a  dose  calculation  module  for  easy  implementation  in  a  radiotherapy  treatment 
planning  system.  We  have  implemented  it  in  an  existing  commercial  RTP  system  for 
conventional  photon/electron  beams  (Ma  et  al  1999a)  and  intensity  modulated  radiotherapy 
(IMRT)  dose  verification  (Ma  et  al  2000a).  MCDOSE  shares  similar  structure  and  geometry 
definition  with  other  EGS4  user  codes,  such  as  DOSXYZ.  MCDOSE  and  DOSXYZ  have 
some  common  features.  They  can  both  simulate  the  transport  of  photons  and  electrons  (and 
positrons)  in  a  3D  rectilinear  phantom  geometry.  The  volume  elements  (voxels)  in  the  phantom 
can  have  uniform  or  variable  dimensions  and  the  material  in  the  voxel  can  be  specified  by  the 
user  or  determined  from  the  electron  density  data  derived  from  the  patient  CT  data.  Source 
models  are  supported  in  both  codes.  But  there  are  other  features  in  MCDOSE,  which  are 
absent  in  DOSXYZ: 

(a)  Advanced  multiple-source  models  (Ma  etal  1997,  Ma  1998,  Deng  etal  2000,  Jiang  et  al 
2000a,  b)  are  used  as  source  input  for  both  photon  and  electron  beams  in  MCDOSE.  The 
multiple-source  models  implemented  in  MCDOSE  are  more  accurate  and  flexible,  and 
easy  for  commissioning. 

(b)  Beam  modifiers  such  as  jaws,  wedges,  blocks,  compensators,  electron  cutouts  and  bolus 
are  included  by  MCDOSE  in  the  patient  simulation.  Both  static  and  dynamic  MLC  fields 
can  be  simulated  for  conventional  and  intensity  modulated  radiotherapy. 

(c)  Several  variance-reduction  techniques  have  been  implemented  in  MCDOSE.  Both  codes 
have  range  rejection.  However,  it  was  implemented  in  different  ways  and  in  different 
simulation  procedures. 

(d)  There  is  an  option  for  selecting  the  geometry  coordinates  from  two  definitions  in 
MCDOSE.  One  (x,  y,  z )  is  the  same  as  that  in  DOSXYZ  for  phantom  study.  The  other 
(x\  /,  z!)  is  for  convenient  patient  dose  calculation  corresponding  to  the  treatment  system. 
The  relationship  between  the  two  coordinates  is  x  —  x'y  y  =  —zf  and  z  =  /- 

(e)  Beamlets  dose  calculation  for  Monte  Carlo  inverse  planning  can  be  performed  by 
MCDOSE  for  both  photon  and  electron  beams. 

(f)  MCDOSE  can  produce  dose  volume  histograms  (DVH)  using  the  patient  contour 
information. 

The  inclusion  of  beam  modifiers  and  the  implementation  of  variance  reduction  techniques 
makes  MCDOSE  more  practical  for  routine  dose  calculation.  Beam  modifiers  play  an 
important  role  in  radiotherapy  treatment.  They  modify  the  beam  shape  and/or  intensity 
distribution  and  therefore  change  the  dose  distribution  in  the  patient.  They  need  to  be 
considered  for  treatment  planning  dose  calculations.  It  is  not  sufficient  to  treat  jaws  and  blocks 
as  apertures  and  to  change  the  profile  using  the  wedge  or  compensator  factors  conventionally. 
The  effect  of  beam  modifiers  has  been  studied  by  Schach  von  Wittenau  etal  (2000).  It  was  found 
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that  the  scattered  electrons  generated  in  jaws  could  increase  the  surface  dose  by  several  per  cent 
depending  on  the  beam  energy  and  field  size.  The  photon  transmission  and  bremsstrahlung 
generated  in  the  modifiers  need  to  be  considered.  To  ensure  simulation  accuracy,  these  beam 
modifiers  should  be  simulated  with  accurate  geometry  when  the  Monte  Carlo  method  is  used 
for  routine  dose  calculation.  Some  codes,  like  BEAM  and  PEREGRINE  (Walling  et  al  1998), 
can  be  used  to  simulate  beam  modifiers.  But  BEAM  can  only  simulate  them  separately  from 
patient  dose  calculation.  One  has  to  simulate  particle  transport  in  the  beam  modifiers  using 
BEAM  first.  After  scoring  the  phase  space  below  the  modifiers,  dose  calculation  in  the  patient 
phantom  can  then  be  performed  (Ma  et  al  1999a).  This  two-step  method  works  well  but  is 
inconvenient  and  can  be  time-consuming.  It  is  not  suitable  for  routine  dose  calculations.  The 
inclusion  of  beam  modifiers  into  the  patient  dose  calculation  simplifies  the  simulation  procedure 
in  terms  of  simulation  steps  and  intermediate  data  storage.  Although  Monte  Carlo  is  an  accurate 
method,  the  CPU  time  needed  to  perform  statistically  meaningful  dose  calculation  is  always  a 
problem  that  prevents  it  from  being  applied  in  routine  dose  calculations.  The  variance  reduction 
techniques  implemented  in  MCDOSE  generally  speed  up  the  dose  calculation  by  a  factor  of 
10-30  (Ma  etal  2000c).  This  makes  the  code  more  practicable  for  routine  application.  The  dose 
calculation  using  MCDOSE  for  a  nine-field  IMRT  plan  (including  pre-  and  post-optimization 
dose  calculation)  can  be  finished  in  1  to  4  h  on  a  450  MHz  Pentium  III  personal  computer 
(Pawlicki  etal  1999,  Ma  etal  2000b). 

Before  MCDOSE  can  be  used  reliably  for  clinical  dose  calculation,  it  must  be  extensively 
validated.  The  accuracy  of  the  dose  calculation  depends  on  the  input  beam  data  and  the 
implementation  of  the  user  code.  Source  modelling  and  beam  commissioning  have  been 
investigated  by  Ma  et  al  (1997),  Ma  (1998),  Deng  et  al  (2000)  and  Jiang  et  al  (2000a,  b). 
Excellent  agreement  in  dose  distributions  in  both  homogeneous  and  heterogeneous  phantoms 
has  been  achieved  between  source  model,  phase  space  data  and  measurements.  In  this  work,  we 
focus  on  the  validation  of  the  user  code  for  beam  modifiers  and  variance  reduction  techniques. 
The  geometry  description  for  beam  modifiers  and  the  implementation  of  variance  reduction 
techniques  may  affect  the  simulation  accuracy  of  the  code. 

The  validation  of  beam  modifiers  was  performed  by  comparing  the  simulation  results  of 
MCDOSE  with  those  of  BEAM/DOSXYZ.  The  dose  distributions  in  specifically  designed 
homogeneous  and  inhomogeneous  phantoms  with  different  source  inputs  calculated  using 
MCDOSE  and  DOSXYZ  are  presented.  Machine  data  for  two  Varian  accelerators,  Clinac 
2100C  and  2300C/D,  were  measured  as  required  by  the  commissioning  procedure. 

2.  Materials  and  methods 

2.7.  Simulation  of  beam  modifiers  in  MCDOSE 

In  this  work,  we  have  implemented  beam  modifier  simulation  in  MCDOSE.  The  combination 
of  the  modifiers,  the  geometry  parameters  and  materials  can  be  input  by  the  user.  Several 
different  modifiers  can  be  used  together.  Each  of  them  occupies  a  slab.  For  photon  beams  our 
multiple-source  model  covers  all  the  fixed  components  above  the  jaws  in  the  accelerator  head 
(Deng  et  al  2000).  What  we  need  to  consider  for  routine  dose  calculations  are  the  jaws,  wedges, 
blocks,  compensators,  MLC  and  bolus.  For  electron  beams  our  multiple- source  model  covers 
all  the  components  down  to  the  lowest  scraper  of  the  applicator  (Jiang  et  al  2000a).  Only  the 
electron  cutout  and  bolus  need  to  be  considered  in  MCDOSE. 

The  x  and  y  jaws  are  modelled  as  two  pairs.  The  inner  surface  of  each  jaw  focuses  to  the 
target.  The  input  parameters  are  x  1 ,  x2,  y  1 ,  y2  and  the  material  for  the  jaws.  Here  x  1 ,  x2 ,  y  1 , 
y2  specify  the  locations  of  the  inner  edges  of  the  jaw  bars.  They  correspond  to  the  field  size 
at  100  cm  SSD.  A  schematic  diagram  for  the  jaws  is  shown  in  figure  1(a). 
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Figure  1.  Schematic  diagram  of  the  beam  modifiers  simulated  in  MCDOSE.  ( a )  X  and  Y  jaws. 
( b )  Photon  beam  block  with  tray,  (c)  Wedge  and  tray.  ( d )  Compensator  and  tray. 


To  simulate  the  treatment  block  and  the  tray,  the  location,  thickness  and  material  of  the 
block  and  the  tray  are  required.  The  location  is  defined  as  the  distance  to  the  isocentre  from  the 
bottom  surface  of  the  block  tray.  The  user  also  needs  to  specify  the  opening  by  the  coordinates 
of  the  vertices  that  are  projected  at  100  cm  SSD.  The  planes  defining  the  inner  surfaces  of  the 
opening  are  all  angled  with  respect  to  the  Z-axis  towards  the  target  as  a  single  focus  point.  The 
block  opening  can  be  of  any  shape.  The  user  needs  to  input  the  points  continuously  around 
the  perimeter  of  the  opening.  The  beam’s  eye  view  and  X-Z  view  of  the  block  and  tray  are 
shown  in  figure  1(b). 

Wedges  can  be  simulated  by  MCDOSE  using  the  geometry  information  about  the 
orientation,  location,  material  and  the  two-dimensional  point  coordinates  to  specify  its  shape 
(see  figure  1(c)).  Full  simulation  of  the  particle  transport  is  performed.  This  approach 
is  similar  to  the  wedge  definition  for  conventional  dose  calculation  algorithms,  except  for 
one  important  point.  For  the  Monte  Carlo  approach,  the  user  enters  the  exact  wedge 
specifications  without  any  modifications  or  adjustments  in  the  parameters  so  that  the  final 
dose  calculation  matches  the  measured  dose  profiles  of  the  wedge  under  consideration.  Thus, 
in  our  Monte  Carlo  approach,  the  physical  transport  of  particles  in  the  beam  modifiers  and  the 
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reproduction  of  the  phase  space  must  be  accurate  to  obtain  good  agreement  with  measured 
results. 

The  compensator  including  the  tray  can  be  simulated  by  MCDOSE.  The  location,  material 
and  the  tray  thickness  are  required.  The  compensator  is  divided  into  different  bins  in  two 
directions.  The  coordinates  and  thickness  for  each  bin  (a  pixel)  are  required.  An  X-Z  view 
of  the  compensator  and  tray  is  shown  in  figure  1  (d).  This  module  can  also  be  used  as  an  MLC 
after  setting  the  tray  material  to  air. 

An  electron  cutout  is  simulated  in  a  way  similar  to  a  photon  block  except  that  it  does  not 
have  a  tray.  The  option  to  have  a  diverging  or  straight  and  parallel  inner  planes  to  define  the 
opening  is  considered  for  different  clinical  applications.  The  location,  material  and  coordinates 
of  each  vertex  to  specify  the  opening  are  required. 

The  bolus  is  simulated  by  adding  an  extra  layer  of  material  to  the  patient’s  geometry 
(outside  the  patient  external  contour)  according  to  its  material  and  thickness. 

2.2.  Simulation  of  beam  modifiers  by  EGS4/BEAM 

In  order  to  validate  the  code  for  beam  modifier  simulation  in  MCDOSE,  EGS4/BEAM  was  used 
to  simulate  block,  wedge,  compensator  and  electron  cutout.  The  simulations  started  with  the 
phase  space  above  the  jaws.  No  variance  reduction  techniques  were  used  in  these  simulations. 
The  cut-off  energy  of  electron  and  photon  were  set  to  be  700  keV  and  10  keV  (total  energy). 
After  the  phase  space  data  were  obtained  below  the  modifiers,  DOSXYZ  was  used  to  carry 
out  the  dose  calculation.  Comparison  between  them  can  prove  that  our  implementation  for  the 
simulation  of  beam  modifier  is  right. 

Two  component  modules,  BLOCK  and  SLABS,  were  used  to  simulate  the  photon  blocks 
and  the  trays.  To  simulate  the  wedge,  the  left  bars  of  two  or  more  pairs  of  JAWS  were  piled 
up  to  match  the  shape.  The  component  modules  MLC  and  SLABS  were  used  to  simulate  the 
compensators  and  trays.  The  component  BLOCK  was  used  to  simulate  an  electron  cutout  with 
the  phase  space  above  the  last  scraper  as  source  input.  Five  subregions  were  used  to  specify 
a  butterfly  shape  for  the  opening  that  was  simulated  for  an  electron  cutout  in  this  paper.  The 
phase  space  data  were  scored  below  the  modifiers.  Then  dose  calculations  were  performed 
using  DOSXYZ  in  a  water  phantom. 

2.3.  Variance  reduction  techniques  and  dose  calculation  details 

The  goal  of  MCDOSE  is  to  perform  quick  and  accurate  dose  calculation  for  radiation 
therapy  treatment  planning.  Several  variance  reduction  techniques  have  been  implemented 
in  MCDOSE  to  speed  up  the  calculation  without  losing  accuracy.  These  techniques  include 
electron  range  and  region  rejection,  photon  interaction  forcing,  particle  splitting,  Russian 
roulette  (Bielajew  and  Rogers  1988)  and  electron  track  repeating  (Kawrakow  etal  1996,  Keall 
and  Hoban  1996b).  To  speed  up  the  simulation  in  beam  modifiers,  the  electron  region  rejection 
technique  was  applied  to  the  electrons  that  make  little  contribution  to  the  dose.  For  example, 
when  we  transport  particles  in  jaws  and  blocks,  we  have  an  option  to  transport  electrons  only 
in  a  margin  around  the  opening.  Variable  global  cut-off  energy  (ECUT)  has  been  applied  to 
different  regions.  The  techniques  of  photon  interaction  forcing,  splitting  and  Russian  roulette 
are  well  known  and  are  implemented  in  some  codes,  such  as  BEAM/EGS4  system.  They  were 
well  described  by  Bielajew  and  Rogers  (1988)  and  Rogers  etal  (1995).  Electron  track  repeating 
is  similar  to  a  technique  called  correlated  sampling  (Ma  and  Nahum  1993),  in  which  particle 
tracks  were  repeated  in  different  locations  or  geometries  to  improve  computation  efficiency.  It 
was  described  in  great  detail  by  Kawrakow  et  al  (1996)  and  Keall  and  Hoban  (1996b).  The  main 
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Figure  2.  Comparisons  of  the  lateral  dose  profiles  at  different  depths  in  a  water  phantom  between 
MCDOSE  and  BEAM/DOSXYZ.  A  45°  wedge  were  simulated  for  a  15  MV  photon  beam,  the  field 
size  was  10  cm  x  10  cm  defined  at  100  cm  SSD. 


idea  is  to  simulate  an  electron  in  a  standard  phantom  such  as  in  water  or  tissue,  recording  all 
the  information  of  the  electron  transport  tracks  including  the  step  length,  direction,  energy 
deposition,  type  of  collision  and  generated  particles  for  each  step,  then  repeat  the  tracks  for 
other  electrons  with  the  same  energy  starting  from  a  different  location  in  the  patient.  The 
initial  direction  for  each  electron  may  be  different.  The  length  of  every  step  in  the  electron 
track  is  inversely  proportional  to  the  material  density  at  that  location.  Stopping  power  ratio  and 
scattering  power  ratio  for  different  materials  need  to  be  considered  to  scale  the  path  lengths 
and  multiple  scattering  angles.  This  technique  saves  the  sampling  time  for  multiple  scattering 
and  collision  for  all  the  repeated  electrons  but  not  for  the  initial  one.  How  many  times  to  repeat 
the  tracks  depends  on  the  problem.  The  number  may  affect  the  final  efficiency.  Our  goal 
was  to  integrate  this  technique  efficiently  with  other  techniques  mentioned  above  to  produce 
a  fast  Monte  Carlo  code.  For  example,  when  photons  go  into  the  phantom,  they  are  split  and 
forced  to  have  same  interactions  at  different  locations.  Many  electrons  are  generated  and  their 
tracks  are  repeated  at  these  locations.  So  time  is  saved  because  we  do  not  perform  initial 
simulation  for  all  the  electrons.  Russian  roulette  is  applied  to  reduce  the  number  of  scattering 
photons. 

To  ensure  the  accuracy  of  the  dose  calculation  with  MCDOSE,  specially  designed 
inhomogeneous  phantoms,  such  as  layered-lung  or  layered-bone  phantoms,  are  used.  The 
dose  distributions  (depth  dose  curves  and  lateral  dose  profiles)  calculated  by  MCDOSE  and 
DOSXYZ  for  different  photon  beams  and  electron  beams  are  compared.  DOSXYZ  was 
selected  for  comparison  with  MCDOSE  because  it  has  already  been  tested  thoroughly  and 
it  is  also  an  EGS4  user  code.  Agreement  between  MCDOSE  and  DOSXYZ  means  that  the 
beam  modifier  implementation  in  MCDOSE  is  correct  and  the  application  of  variance  reduction 
techniques  does  not  affect  the  accuracy.  Since  DOSXYZ  has  been  well  benchmarked  by  other 
investigators  and  our  previous  work  (Rogers  et  al  1995,  Kapur  et  al  1998,  Ma  1998,  Zhang 
et  al  1999,  Ma  et  al  1999a,  2000a),  it  is  not  necessary  to  compare  all  the  MCDOSE  results 
with  measurements.  In  this  work,  both  monoenergetic  photon  beams  and  realistic  clinical 
beams  of  Varian  2300C/D  are  studied.  The  goal  for  use  of  monoenergetic  photon  beams  is  to 
make  sure  that  the  code  works  well  for  monoenergetic  beams  (therefore  it  should  also  work 
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Figure  3.  Dose  distributions  in  a  water  phantom  calculated  by  MCDOSE  (thin  curve)  and 
BEAM/DOSXYZ  (thick  curve)  for  a  10  cm  x  10  cm  field  of  a  15  MV  photon  beam.  The  beam  was 
modified  by  a  hexagonal  block  and  a  60°  wedge,  (a)  The  isodose  distribution  in  the  X-Z  plane 
through  the  central  axis  (beams  come  from  the  bottom).  ( b )  The  isodose  distribution  in  the  X-Y 
plane  at  a  depth  of  3.5  cm. 


for  polyenergetic  beams).  This  test  was  done  thoroughly  in  the  energy  range  of  interest  for 
dosimetry  with  an  interval  of  1  MeV. 

3.  Results 

3.1.  Validation  of  beam  modifiers 

3.1.1.  A  45°  wedge.  Figure  2  shows  the  lateral  dose  profiles  at  different  depths  in  a  water 
phantom  calculated  using  MCDOSE  and  BEAM/DOSXYZ  for  a  45°  wedge.  Lead  was  used 
as  the  wedge  material.  The  field  size  was  set  to  10  cm  x  10  cm  at  100  cm  SSD.  The  voxel 
size  in  the  water  phantom  was  1.0  cm  x  0.3  cm  x  0.5  cm.  A  15  MV  photon  beam  from  a 
Varian  Clinac  2300C/D  machine  was  used  as  source  input  with  phase  space  scored  above  the 
jaws.  The  comparison  calculations  were  performed  under  the  same  conditions.  The  agreement 
between  the  results  of  MCDOSE  and  BEAM/DOSXYZ  is  within  1%  of  the  maximum  dose. 
The  CPU  time  comparison  between  MCDOSE  and  BEAM/DOSXYZ  is  shown  in  table  1 . 
MCDOSE  is  about  16  times  faster  than  BEAM/DOSXYZ  for  this  simulation. 

3.1.2.  A  hexagonal  block  and  a  60°  wedge.  Figure  3  shows  the  dose  distributions  calculated  by 
MCDOSE  and  BEAM/DOSXYZ  in  a  water  phantom  with  a  hexagonal  block  and  a  60°  wedge 
for  a  10  cm  x  10  cm  field  of  a  1 5  MV  photon  beam.  The  block  was  made  of  Cerrobend  and  its 
thickness  was  7.7  cm.  The  block  tray  was  made  of  PMMA  and  its  thickness  was  0.7  cm.  The 
wedge  was  made  of  lead.  The  voxel  size  in  the  water  phantom  was  0.5  cm  x  0.5  cm  x  0.5  cm. 
Figure  3(a)  shows  the  isodose  distribution  in  the  X-Z  plane  through  the  central  axis  and 
figure  3(b)  shows  the  isodose  distribution  in  the  X-Y  plane  at  a  depth  of  3.5  cm.  The  difference 
between  the  results  of  MCDOSE  and  BEAM/DOSXYZ  is  within  1%  of  the  maximum  dose. 
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Figure  4.  (a)  The  isodose  distributions  in  XT  view  at  a  depth  of  2.5  cm  in  a  water  phantom  for 
a  modified  10  cm  x  10  cm  15  MV  photon  beam  using  a  specially  designed  compensator.  The 
thick  curve  is  the  simulation  result  of  BEAM/DOSXYZ.  The  thin  curve  is  the  result  of  MCDOSE. 
(b)  The  lateral  dose  profile  along  the  dashed  line  in  (a).  The  curve  with  symbols  in  ( b )  is  the  result 
of  MCDOSE. 


Table  1 .  The  CPU  time  (hours)  required  by  MCDOSE  and  BEAM/DOSXYZ  for  the  beam  modifier 
simulation  investigated  in  this  work  on  a  Pentium  III  450  MHz  personal  computer.  The  statistical 
uncertainties  of  the  maximum  doses  are  1%. 


Modifier 

45°  wedge 

60°  wedge/block 

E-shaped 

compensator 

Butterfly-shaped 
electron  cutout 

BEAM  +  DOSXYZ 

MCDOSE 

Ratio 

31.5  +  6.1 

2.3 

16.3 

247.5  +  7.9 

13.0 

19.4 

10.3+10.9 

5.3 

4.0 

37.9  +  1.4 

2.5 

15.7 

The  CPU  time  needed  by  MCDOSE  and  BEAM/DOSXYZ  is  given  in  table  1.  MCDOSE  is 
about  19  times  faster  than  BEAM/DOSXYZ  for  this  simulation. 
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Figure  5.  ( a )  is  the  isodose  distributions  in  XY  view  at  a  depth  of  2.25  cm  in  a  water  phantom  for 
a  modified  10  cm  x  10  cm  12  MeV  electron  beam  using  an  electron  cutout  with  a  butterfly-shaped 
opening.  The  thick  curve  is  the  simulation  result  of  BEAM/DOSXYZ.  The  thin  curve  is  the  result 
of  MCDOSE.  (fr)  The  lateral  dose  profile  along  the  dashed  line  in  (a).  The  curve  with  symbols  in 
( b )  is  the  result  of  MCDOSE. 


3.1.3.  Compensator.  Figure  4  shows  the  E-shaped  compensator  isodose  distributions  (a)  and 
lateral  dose  profiles  ( b )  at  a  depth  of  2.5  cm  in  the  water  phantom.  They  were  calculated  using 
MCDOSE  and  BEAM/DOSXYZ  respectively.  The  compensator  was  made  of  copper.  The 
thickness  in  the  area  of  the  compensator  forming  the  ‘E’  shape  was  0.5  cm  and  the  thickness  in 
the  remaining  area  was  2  cm.  The  voxel  size  in  the  water  phantom  was  0.4  cm  x  0.4  cm  x  1 .0  cm. 
Agreement  (within  1%)  can  be  observed  from  the  comparison.  In  the  penumbra  regions,  the 
separation  between  the  results  of  MCDOSE  and  BEAM/DOSXYZ  is  within  1  mm.  The  CPU 
time  needed  by  MCDOSE  and  BEAM/DOSXYZ  is  also  shown  in  table  1 .  MCDOSE  is  about 
four  times  faster  than  BEAM/DOSXYZ  for  this  simulation. 

3.1.4.  Electron  cutout.  Figure  5  shows  the  isodose  distributions  {a)  and  lateral  dose 
profiles  ( b )  in  a  water  phantom  calculated  using  BEAM/DOSXYZ  and  MCDOSE,  for  a  1 2  MeV 
realistic  electron  beam  from  a  Varian  21 00C  machine  with  a  butterfly-shaped  cutout.  The  cutout 
was  made  of  Cerrobend  and  its  thickness  was  1.4  cm.  The  voxel  size  in  the  water  phantom 
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Figure  6.  Dose  distribution  in  a  tissue-lung-tissue  phantom  calculated  by  DOSXYZ  and  MCDOSE 
for  6  MeV,  15  MV  and  2  McV  photon  beams  with  field  size  of  6  cm  x  6  cm  defined  at  100  cm 
SSD.  The  material  in  the  volume  of  —5.0  cm  to  5.0  cm  along  the  X  direction,  —2.0  cm  to  2.0  cm 
along  the  Y  direction  and  7.0  cm  to  12.5  cm  along  the  Z  direction  is  lung,  (a)  The  depth  dose 
curves  along  the  central  axis  for  the  beams.  ( b )  The  lateral  dose  profiles  along  the  X-axis  at  depths 
of  3.3  cm,  9.0  cm  and  15.0  cm  for  the  15  MV  photon  beam. 


was  0.5  cm  x  0.5  cm  x  0.3  cm.  The  phase  space  of  an  electron  beam  above  the  last  scraper 
of  a  10  cm  x  10  cm  applicator  was  used  as  the  source  input.  Good  agreement  (within  1%) 
has  been  achieved  between  the  results  of  MCDOSE  and  BEAM/DOSXYZ.  In  the  penumbra 
regions,  the  separation  between  MCDOSE  and  BEAM/DOSXYZ  is  within  1  mm.  The  CPU 
time  needed  by  MCDOSE  and  BEAM/DOSXYZ  is  shown  in  table  1.  MCDOSE  is  about  16 
times  faster  than  BEAM/DOSXYZ  for  this  simulation. 

The  electron  beam  cutout  factors  calculated  using  MCDOSE  are  compared  with  the 
published  data  in  table  2.  The  published  data  were  calculated  using  BEAM/DOSXYZ  with  full 
phase  space  input  and  measured  by  Kapur  et  al  (1998).  The  difference  between  the  MCDOSE 
data  and  those  calculated  with  full  phase  space  by  DOSXYZ  is  within  1 .0%.  The  differences 
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Figure  7.  Dose  distributions  in  a  tissue-bone-tissue  phantom  calculated  by  DOSXYZ  and 
MCDOSE  for  6  MeV,  15  MV  and  2  MeV  photon  beams  with  field  size  of  6  cm  x  6  cm  defined  at 
100  cm  SSD.  The  material  in  the  volume  of  —5.0  cm  to  5.0  cm  along  the  X  direction,  —2.0  cm  to 
2.0  cm  along  the  Y  direction  and  7.0  cm  to  10.0  cm  along  the  Z  direction  is  bone.  ( a )  The  depth 
dose  curve  along  the  central  axis  for  the  beams.  ( b )  The  lateral  dose  profiles  along  the  Y  axis  at 
depths  of  3.3  cm,  8.0  cm  and  1 1.0  cm  for  the  15  MV  photon  beam. 

in  the  cutout  factors  between  MCDOSE  data  and  measurement  is  less  than  2%  except  for  the 
cases  of  20  MeV  electron  beams  with  3  cm  x  3  cm  and  4  cm  x  4  cm  inserts. 


3.2.  Validation  of  dose  distribution  in  heterogeneous  phantoms 

To  test  the  accuracy  of  the  dose  distribution  calculated  using  MCDOSE  with  the  variance 
reduction  techniques,  comparisons  are  made  with  the  results  of  DOSXYZ.  Each  example  has 
been  calculated  with  a  source-surface  distance  of  100  cm  and  field  size  of  6  cm  x  6  cm. 
The  global  electron  cut-off  energy  (ECUT)  was  set  to  700  keV  and  photon  cut-off  energy 
(PCUT)  was  set  to  10  keV  in  MCDOSE  and  DOSXYZ.  The  maximum  fractional  energy  loss 
per  electron  step  (ESTEPE)  was  limited  to  4%  for  both  codes.  The  material  compositions  and 
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Figure  8.  Dose  distributions  in  a  tissue-lung-tissue  phantom  calculated  by  DOSXYZ  and 
MCDOSE  for  a  15  MeV  electron  beam  with  field  size  6  cm  x  6  cm.  The  material  in  the  volume 
of  —5.0  cm  to  5.0  cm  along  the  X  direction,  —2.0  cm  to  2.0  cm  along  Y  direction  and  4.2  cm  to 
6.0  cm  along  the  Z  direction  is  lung,  (a)  The  depth  dose  curve  along  the  central  axis.  ( b )  The 
lateral  dose  profiles  along  the  K-axis  at  depths  of  3.0  cm,  5.0  cm  and  7.0  cm. 


densities  are  taken  from  ICRU  (1992).  Monoenergetic  point  sources  and  phase  space  files  for 
a  15  MV  photon  beam  of  a  Varian  Clinac  2300C/D  linear  accelerator  were  used.  Since  dose 
per  incident  fluence  is  calculated,  absolute  dose  comparisons  are  made  between  MCDOSE  and 
DOSXYZ. 

3.2.1.  Photon  beams.  Figure  6(a)  shows  the  depth  dose  curves  along  the  central  axis  for 
6  MeV  and  2  MeV  monoenergetic  photon  beams  calculated  by  DOSXYZ  and  MCDOSE  in  a 
tissue  phantom  with  3D  lung  inhomogeneity.  At  depths  from  7  cm  to  12.5  cm,  the  phantom 
contains  a  10  cm  x  4  cm  slab  of  lung  in  the  centre.  The  depth  dose  curves  for  a  realistic 
15  MV  photon  beam  from  a  Varian  Clinac  2300C/D  accelerator  calculated  by  DOSXYZ  and 
MCDOSE  in  the  same  phantom  are  also  shown  in  figure  6(a)  and  the  lateral  dose  profiles  at 
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Figure  9.  Dose  distribution  in  tissue-bone-tissue  phantoms  calculated  by  DOSXYZ  and  MCDOSE 
for  15  MeV  and  6  MeV  electron  beams  with  field  size  6  cm  x  6  cm.  For  the  15  MeV  electron 
beam,  the  material  in  the  volume  of  —5.0  cm  to  5.0  cm  along  the  X  direction,  —2.0  cm  to  2.0  cm 
along  the  Y  direction  and  2.7  cm  to  4.5  cm  along  the  Z  direction  is  bone.  For  the  6  MeV  electron 
beam,  the  material  in  the  volume  of  -5.0  cm  to  5.0  cm  along  the  X  direction,  -2.0  cm  to  2.0  cm 
along  the  Y  direction  and  0.6  cm  to  1 .6  cm  along  the  Z  direction  is  bone,  (a)  The  depth  dose  curves 
along  the  central  axis  for  two  beams,  (b)  The  lateral  dose  profiles  along  the  Y  axis  at  depths  of 
3.4  cm  and  5.0  cm  for  the  15  MeV  electron  beam. 

different  depths  are  shown  in  figure  6(b).  The  phase  space  for  a  15  MV  photon  beam  was 
generated  by  BEAM  and  contains  over  5  million  particles.  The  lcr  statistical  uncertainty  on 
the  dose  values  for  these  curves  is  less  than  0.5%. 

Dose  calculations  were  also  performed  for  photon  beams  in  an  inhomogeneous  phantom 
containing  tissue  and  bone.  A  1 0  cm  x  4  cm  slab  of  bone  was  placed  between  7  cm  and  10  cm 
depth  in  the  centre  of  the  tissue  phantom.  Figure  1(a)  shows  the  depth  dose  curves  along  the 
central  axis  for  different  photon  beams.  Figure  1(b)  shows  the  lateral  dose  profiles  at  different 
depths  for  a  15  MV  photon  beam. 

Figures  6  and  7  show  good  agreement  in  dose  distributions  between  MCDOSE  and 
DOSXYZ  for  different  photon  beams.  The  absolute  doses  calculated  by  these  two  codes 
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Table  2.  The  electron  cutout  factors  for  various  square  inserts  in  10  x  10  cm2  applicator  for  6,  12 
and  20  MeV  beams  calculated  by  MCDOSE  and  compared  with  the  full  phase  space  calculation 
by  DOSXYZ  and  measurement  of  Kapur  et  al  (1998).  The  values  in  the  parentheses  indicate  the 
per  cent  difference  between  the  MCDOSE  calculated  data  and  those  measured  or  calculated  with 
full  phase  space  by  DOSXYZ. 


Energy  (MeV) 

Insert  (cm  x  cm) 

MCDOSE 

Full  phase  space 

Measurement 

6 

3x3 

0.922 

0.923  (-0.1%) 

0.927  (-0.5%) 

4x4 

0.983 

0.982  (0.1%) 

0.988  (-0.5%) 

8x8 

1.009 

1.005  (0.4%) 

1.003  (0.6%) 

12 

3  x  3 

0.930 

0.930  (0.0%) 

0.928  (0.2%) 

4x4 

0.965 

0.956  (0.9%) 

0.963  (0.2%) 

8x8 

0.999 

1.002  (-0.3%) 

0.991  (0.8%) 

20 

3x3 

0.964 

0.968  (-0.4%) 

0.993  (-2.9%) 

4x4 

0.989 

0.993  (-0.4%) 

1.011  (-2.2%) 

8x8 

0.999 

0.993  (0.6%) 

1.004  (-0.5%) 

in  tissue,  lung,  bone,  even  at  the  interface  between  tissue  and  bone  and  tissue  and  lung  agree  to 
within  1  %  of  the  maximum  dose.  Further  comparisons  of  dose  calculations  for  monoenergetic 
photon  beams  from  1  MeV  to  20  MeV,  and  also  for  6  MV  and  1 0  MV  photon  beams  were 
performed.  Similar  agreement  between  MCDOSE  and  DOSXYZ  was  achieved.  In  most  cases, 
the  differences  were  within  the  lcr  statistical  uncertainty  of  1%. 

3.2.2.  Electron  beams.  We  also  performed  dose  calculations  for  monoenergetic  point  source 
electron  beams  with  field  size  6  cm  x  6  cm.  The  depth  dose  curves  of  a  15  MeV  electron 
beam  calculated  by  MCDOSE  and  DOSXYZ  in  a  tissue-lung-tissue  phantom  are  shown  in 
figure  8(a).  The  lung  material  was  at  a  depth  from  4  to  6  cm.  Figure  8 (b)  shows  the  lateral  dose 
profiles  at  different  depths  in  the  phantom.  At  a  depth  from  4.2  to  6.0  cm,  the  lung  material 
was  from  —5  cm  to  5  cm  in  the  X  direction  and  from  —2  cm  to  2  cm  in  the  Y  direction.  The 
remaining  volume  was  tissue.  The  difference  between  the  simulation  results  of  MCDOSE  and 
DOSXYZ  is  within  1%  of  the  maximum  dose. 

The  depth  dose  curves  and  lateral  dose  profiles  at  different  depths  for  a  15  MeV  electron 
point  source  beam  in  a  tissue-bone-tissue  phantom  are  shown  in  figures  9(a)  and  (b).  At  a 
depth  from  2.7  to  4.5  cm,  there  was  a  10  cm  x  4  cm  slab  of  bone  in  the  centre  surrounded  by 
tissue.  The  depth  dose  curve  for  a  6  MeV  electron  beam  in  a  tissue-bone-tissue  phantom  is 
also  shown  in  figure  9(a).  The  bone  material  was  at  a  depth  from  0.6  to  1 .6  cm.  The  difference 
between  the  simulation  results  of  MCDOSE  and  DOSXYZ  is  within  1  %  of  the  maximum  dose. 

3.2.3.  CPU  time  comparison.  The  CPU  times  required  by  MCDOSE  and  DOSXYZ  for  dose 
calculations  shown  in  figures  6  and  9  are  listed  in  table  3  for  different  photon  and  electron  beams. 
The  calculations  were  performed  on  a  Pentium  III  450  MHz  PC.  The  statistical  uncertainty 
of  the  maximum  dose  is  1%.  The  ratios  in  the  table  are  the  speed-up  factors  of  MCDOSE 
compared  with  DOSXYZ  for  each  beam.  In  general  MCDOSE  is  20-50  times  faster  than 
DOSXYZ  for  photon  beams  and  5-10  times  faster  than  DOSXYZ  for  electron  beams. 


4.  Summary 

A  new  Monte  Carlo  dose  calculation  tool  has  been  developed  and  validated  for  clinical  dose 
calculation.  The  main  features  of  MCDOSE  include  electron  and  photon  beam  reconstruction 
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Table  3.  The  CPU  time  (hours)  required  by  MCDOSE  and  DOSXYZ  for  the  dose  calculations 
shown  in  figures  6  and  9  on  a  Pentium  III  450  MHz  personal  computer.  The  statistical  uncertainties 
of  the  maximum  doses  are  1%. 


Beam 

2  MeV  photon 

6  MeV  photon 

15  MV  photon 

6  MeV  electron 

15  MeV  electron 

DOSXYZ 

9.44 

2.95 

4.86 

0.78 

0.94 

MCDOSE 

0.25 

0.074 

0.098 

0.08 

0.17 

Ratio 

37.7 

39.8 

49.6 

9.8 

5.6 

by  source  models,  simulation  of  beam  modifiers  together  with  the  patient  geometry  for  dose 
calculation  and  application  of  variance  reduction  techniques. 

Comparisons  of  dose  distribution  with  wedges,  blocks,  compensator  and  electron  cutout 
between  MCDOSE  and  BEAM/DOSXYZ  demonstrated  that  dose  calculation  with  modifiers 
in  MCDOSE  is  accurate,  convenient  and  efficient.  To  compare  MCDOSE  and  DOSXYZ  dose 
distributions,  3D  heterogeneous  phantoms  containing  lung  and  bone  were  used.  Excellent 
agreement  between  the  MCDOSE  and  DOSXYZ  results  has  been  obtained  for  both  photon 
and  electron  beams  of  different  energies.  Comparing  with  DOSXYZ,  the  MCDOSE  code  runs 
generally  20-50  times  faster  for  the  photon  beams  and  5-10  times  faster  for  electron  beams 
investigated  in  this  work  without  beam  modifiers. 
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